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Photograph taken at C-I-L Central Research Laboratory, McMasterville, Que. 


Plastics Detective 


This photograph shows the High Shear Viscometer, designed and built at our 


Central Research Laboratory, used to investigate factors involved in complex 
flow behaviour of molten thermo-plastics during commercial processing. 
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Information on this newly developed pressure-driven capillary Viscometer ay 


has been made available to industry by C-I-L, and many industrial laboratories Man 
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in the United States, England, India and Australia have built identica! machines. . a 
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Such C-I-L research contributes to higher standards of performance in the 
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plastics industry, to new and better plastics products and—through better technical ta 
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service—to greater success in the application of these new materials to modern living. Street, ( 
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An Analysis of Spray Evaporation 
in a High-Temperature Environment 


I. Radiant Heat Transfer to Clouds of Droplets and Particles! 


T. W. HOFFMAN? 


This series of papers presents an analysis of the 
rate of evaporation of a spray which is injected into 
a cylindrical column, the walls of which are main- 
tained at a high temperature. The method is an exten- 
sion of the Hottel-Cohen analysis of the radiant heat 
transfer to an absorbing gas of non-uniform tempera- 
ture’! and permits the calculation of the magnitude 
and distribution of the heat flux at the wall. 


The present paper discusses the difficulties in- 
herent in the calculation of the absorption of thermal 
radiation by clouds of particles suspended in an ab- 
sorbing gas medium. Approximations to the intracta- 
ble exact solution are made to permit the calculation 
of the absorption strength of spray from operating 
conditions and the properties of the atomized liquid. 


The various reception factors for the cylindrical 
system have been accurately determined as a fune- 
tion of the combined absorption strength of spray and 
medium. The use of these factors in the present sys- 
tem is discussed. 


hen a cloud of droplets flowing through a cylindrical 

enclosure is subjected to thermal r: adiation, ev aporation 
will occur over at least part of the column. ‘The vapors thus 
generated are generally athermanous, i.e., they will absorb 
thermal radiation in certain wave- length regions of the infra-red 
spectrum. As a result, an additional heat load, other than the 
sensible heat increase of the gas, will be maintained in the 
column, causing severe variations in the temperature of the gas 
from the top — or inlet of the column — to the bottom. 


I his paper deals with the problem of estimating the radiation 
from the heated walls of a cylindrical enclosure to a cloud of 
droplets of particles contained in an absorbing medium, the 
temperature of which varies from top to bottom. In addition, 
the temperature of the walls of the column may vary over the 
column height, depending on the process requirements. Further 
more, the cloud will be assumed to consist ot droplets and /o1 
particles which are in the size range from LO to 200 microns 
Uhese droplets are of a size comparable to the wave lengths 
under consideration, since most of the infra-red spectrum 1s 
contained in the I-to-15 micron wave length region. Moreover, 
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and W. H. GAUVIN: 


Dans cette série d’articles on présente une analyse 
du taux d’évaporation d°un jet pulvérisé injecté dans 
une colonne cylindrique dont les parois sont mainte- 
nues a une température élevée. La méthode décrite 
est une extension de l’analyse de Hottel et Cohen pour 
le transport de la chaleur par rayonnement 4 un gaz 
absorbant, de température non-uniforme. Elle permet 
de calculer lintensité et la distribution du taux de 
transport de chaleur a la paroi. 


Les auteurs exposent les difficultés inhérentes au 
calcul de l'absorption d’énergie rayonnante par des 
suspensions de particules dans un milieu gazeux ab- 
sorbant. Ils ont recours & des approximations des 
solutions exactes, difficiles 4 résoudre, pour caleuler 
le pouvoir absorbant d'un jet pulvérisé 4 partir des 
conditions expérimentales et des propriétés du liquide 
vaporisé, 

On a déterminé de fagon précise les différents 
facteurs de réception pour le systéme cylindrique, en 
fonction du pouvoir absorbant combiné du jet et du 
milieu. On traite de Vutilisation de ces facteurs dans 
le présent systéme. 


these particles may be only partially opaque to the radiation 


striking them. 


l. Gas Radiation 


In a recent paper‘, Horttel and Cohen presented a method 
of calculating the radiant heat exchange in an enclosure containing 
an absorbing medium. In this method, v arying space tempera- 
tures are taken into account. The technique involves dividing 
the enclosure into a number of surface and gas zones, each of 
which can be considered isothermal. ‘Their size and number 
depend upon the accuracy desired. Reception factors corre- 
sponding to Fy, for an on medium in Hottel’s nomencla- 
ture were computed for a square-cube system for surface- 
surface, gas-surtace, gas-gas thermal radiation interchange. The 
determinant technique, presented earlier by Hottel 
of problems in a grey enclosure containing a non-absorbing 
By means of this tech- 


tor solution 


medium, was extended to this system 
nique, total interchange factors corresponding to A 2 can 
be calculated trom the reception factors to allow tor the reflective 
tlux from all surtace zones. Heat and material balances are then 
made on all surtace and gas zones. lo perform this balance on 
a gas zone, for example, one must add together the net radiant 
heat received trom all other gas zones and surtace zones, the 
net convective conductive heat from adjacent gas and surface 
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zones, and the enthalpy increase or decrease in the zone. For 
any such system, the unknowns in the resulting equations are 
the temperatures of the gas zones, and the number of unknown 
temperatures is equal to the number of equations. The solution 
of the simultaneous equations allows a prediction of the temper- 
ature of the gas zones and of the heat flux distribution in the 
enclosure. 

The key parameter in applying this method of solution to 
the radiant heat interchange problem is the absorption coefficient 
k of the gas, which is a measure of the absorption strength of 
the absorbing medium. For a grey gas, & is constant over the 
full wave-length spectrum and the expression relating the 
emissivity or absorptivity of a gas to the absorption coefhicient is: 


geqeat-s? Branarkiar inane a aia w aiehion (1) 
where / is the thickness of the gas. 
However, absorption and emission of thermal radiation by 
a real gas occur only over several well-defined wave- -length 
bands, i.e., only over part of the wave-length spectrum. Hottel ® 
has thus proposed that the actual emissivity or absorptivity 
versus path length curve can be made up of the sum of a number 
of grey gas absorption coefficients, as in Equation (1), but each 
will have a specific weighting factor, namely 


€= Yai (1 ~ chil) where > 4: NE ii xler gel ee (2) 
‘ é 


Proceeding one step further, Hottel suggests that a real gas 
can be replaced by a grey and a clear gas. The grey gas, with 
absorption coefficient k, can be thought of as occupying the 
fraction a of the spectrum and the clear gas (k = 0) occupying 
the fraction 1 — a. In essence, this scheme is one of fitting an 
emissivity versus path length curve to an exponential equation 
of the form: 


€=a(1 — eb) sien ene ions (3) 


It must be remembered that a and & depend on the conditions 
and dimensions of the system under consideration. 


Hottel and Egbert’ » have shown that, for any enclosure in 
which the surface emissivity is fairly high, only the first surface 
reflection of the radiation is of any importance, for the radiation 
undergoes double attenuation on each reflection first by the 
gas and then by the surface. Therefore the emissivity should 

be keyed to Equation (3) at zero path length (/ = 0), at one 
mean beam length L of the enclosure at w hich € = €,, and at 
twice this length (€ = €,). Then a and kL have been shown 
to be given by the relations: 


a = (€,)*/(2e, Gerd. v6 > ~« &#) 
kL = In [€z/( €2z €z)} vies bate (5) 
kZ = (1/W) In [ex /( ee, — €z)). ‘ (6) 


Where Z is the characteristic length of the system and W is 
the factor by which this characteristic length must be multiplied 
to give the mean beam length. Appendix C presents a method 
of finding VY for any cylindrical enclosure containing any given 
concentration of absorbing medium. 

Hottel and Cohen“ present a number of plots of reception 
factors for surface-to-surface, surface to-gas interchange in a 
black enclosure, based on a square and cube system. These 
factors are plotted as a function of kb, where b is the side of a 
cube or square. Total interchange factors can be built up from 
these plots. 

When using these plots and, in fact, when employing this 
method, the absorption coefficient k and the weighting factor a 
are assumed to be constant throughout the enclosure. ‘These 
authors present a detailed diacussion concerning this necessary 
simplification. Their conclusions are that although moderate 
variations may occur in these parameters, the error is minimized, 
for one is sull allowing for the major temperature variation. 
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The factors & and a affect only the coefficients (i.e., the inter- 
change factors); the emissive power (o7;)* is allowed to vary 
at will. 


Thus the problem of solving for temperature and heat flux 
distribution in an enclosure containing a temperature-varying 
absorbing medium can be handled with considerable accuracy, 
However, if one introduces a cloud of particles into this medium, 
the problem becomes more complicated, for now complex surface 
and gas emission or absorption arises in each gas zone. 


2. Radiation to Clouds of Particles 


For the present, the following discussion is restricted to the 
case of radiation to a cloud of particles suspended in a non- 
absorbing medium. 

When thermal radiation traverses a cloud of particles, part 
of it may be absorbed by these particles and part of it may be 
reflected or scattered. The magnitude and direction of the 
scattered portion and the amount “absorbed may be found from 
the equations presented by Mie.* The exact Mie solution 
applies to electro-magnetic waves at any frequency striking a 
particle of any size. 


Both scattering and absorption remove energy from the 
radiant beam in the direction of propagation of the beam.** 
This extinction or attenuation along the direction of travel is 
usually expressed by the Beer-Lambert Law: 


ig ce ae (7) 


where | is the intensity of the radiation after traversing a distance 
/, and J, is the intensity at / = 0. The extinction coefficient 
is computed from the expression 


co 


7 = RR 60 ec wee aan’ (8) 


° 


for the case of particles of different radii. In this expression 
N(r) is the number of particles per unit volume with radii in 
the interval dr. Q... is termed the efficiency factor for extinction 
in the forward direction; it is the ratio of the extinction cross- 
section, S, to the geometric cross-section mr?. It should be 
remembered that Equation (8) is valid only for a relatively 
thin cloud — i.e., for situations where the scattered intensity 
is directly proportional to the number of particles. This implies 
that each particle is essentially exposed to the original radiation. 
With large particle concentrations or relatively great thicknesses, 
multiple scattering may occur. Then this simple proportionality 
does not hold, for now particles are exposed not only to a fraction 
of the original beam, because of extinction, but also to the 
scattered radiation from other particles. This situation arises 
for y/ > 0.1, which is usually the case here. The discussion 
which follows is still valid for the scattering of a single particle, 
but the analysis of the overall picture becomes very complex, 
especially so in those cases where the cloud is restricted to 
finite dimensional limits. At still greater concentrations, i.¢., 
when inter-particle distances are less than 37, phase interactions 
and interference phenomena result, 

The Mie Theory predicts the extinction cross-section, S$, 
for transparent and absorbing spheres. For non-absorbing spheres 


“> (Gg Og) /(20 42) o ccc ccicasn (9) 


n= 1 


where a, and b, are complex functions of the parameter x = 2ar/X 
and the refractive index m, When there is appreciable absorption 
within the sphere, this may be taken into account by employing 


*An excellent presentation of the Mie Theory in condensed form is given 
by V. K. LaMer™, 


°°*Scattered radiation” refers to any radiation which does not proceed 
along the line of original propagation. In this sense, it includes the 
reflected, the refracted and the diffracted part of the beam. 
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a complex refractive index m = n — in’, where n’ is the absorp- 
tion index and is usually a strong function of wave-length. For 
absorbing spheres, Equation (9) gives the scattering cross-section 
only, and the total radiation abstracted from the incident beam 
(scattered plus absorbed) is found from: 


Sear = MS he > (-1)" (an + ba pene 
27 n=1 


where Re stands for the real part of the expression in the brackets. 

Thus, for the problem of thermal radiation to a cloud of 
droplets or particles, the difficulties in applying the Mie solution 
become immediately apparent. These difficulties arise because 
the extinction cross-section, and thus the exponent ‘y in Equation 
(7) depends on 7, \ and m. Here the full spectrum of thermal 
radiation wave-lengths, the distribution of particle size in the 
cloud, and the variable refractive index must be considered. 
Moreover, in practice thick dense clouds are usually encountered 
and multiple scattering results. 

Other difficulties arise when an attempt is made to apply 
the Mie solutions to the problem at hand. Thus, a simplification 
is necessary for engineering applications which is fortunately 
possible under certain conditions. When the parameter x=2m7r/X 
becomes large, the extinction efficiency Q... tends to a limiting 
value of two, that is, the extinction cross-section becomes 277°. * 
Geometrical optics predicts that the covering capacity of a 
spherical particle should be equal to the cross-sectional area 
of the particle, rr. This seeming anomaly has been explained 
by Sinclair and Lothian and Chappel“, who pointed out that 
the Mie Theory accounts for diffraction as well as absorption 
and reflection of the radiant beam by the particle. On the basis 
of Babinet’s principle of diffraction of light, it is shown that a 
circular opaque disc diffracts an amount of light equal to that 
which actually strikes the disc. Furthermore, Gumprecht and 
Sliepcevich“? have shown that this diffracted radiation is 
concentrated within a relatively small cone angle in the forward 
direction (ca. 2° for values of x = 100 and m = 1.33). This 
angle becomes larger as x decreases, but at all values of x that 
are of interest, and certainly for the enclosure under considera- 
tion, this part of the beam need not be considered scattered. 
Thus the main concern is that part of the beam which is reflected 
and/or. absorbed. 

If the average absorptivity of the particle a, is close to one, 
that is, the particles behave very nearly like black bodies, then 
the reflected part of the radiation is inconsequential and for a 
cloud of particles, Equation (7) becomes: 


Fw Fe ces cece s Git) 


where r, is the mean surface radius calculated from: 


= (don ri2)/yon ee . (12) 
é 4 


For this case, Equation (11) is exact; however if the particles 
are only partially absorbing, then this expression (with a, 
removed) is restricted to the direction of propagation only. 
Equation (11) has been derived previously by Haslam and 
Hottel from geometrical optics considerations and was used 
by Johnstone et al“ in their studies on heat transfer to clouds 
of particles. 

he absorptivity of a droplet depends on its material compo- 
sition, the wave length of radiation striking it, and upon its 
size (more exactly its volume). As mentioned previously, it 
is possible to compute the exact absorption within the particle 
by employing the Mie equations, It is interesting to note that 
for particles approaching black conditions, a minimum size is 


*The scattering cross-section as calculated from the Mie Theory is not 
2Wr* for all x. It can be more or less than 27r*® depending upon the 
interference effects of the reflected, refracted and diffracted radiation. 
The laws of geometrical optics are asymptotic laws for electromagnetic 

waves, Geometrical or rayoptic methods may be used when the particles 
become large in comparison to the wave-length of the radiation striking 
the particle, 
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indicated by these equations. Van de Hulst“® points out that 
two conditions are necessary: 


(1) that a negligible fraction is reflected; 
(2) that the sphere is not transparent. 


Condition (1) implies that |—1! << 1 so that m and n’ must 
be small. If is not close to 1, reflection at the surface according 
to Fresnel’s formula is: 

1] sin? (6 — p) 
2} sin? (6 + W) 


tan? (@ — p) 


Ro => Oe 
tan? (6 + yp) 


(13) 


where Rg is the fraction reflected and from the refractive law 
(see Figure 1(a)): 


SO MGs. oi cack ce ceenet (14) 


Condition (2) implies that n’ x >> 1 and since condition (1) 
indicated that 7’ is small, then x = 247r/X must be large (ca. 100, 
or the radius of the particle must be larger than about fifteen 
times the wave-length). 

To circumvent the complex Mie solution, geometrical optical 
methods can be applied to this absorption problem. Nusselt“ 
has shown that the absorptivity of a sphere containing an absorb- 
ing medium with an average absorption coefficient K can be 
calculated from: 


= 1 + [e-*®/Kr] — [1 — e-®®]/2(Kr)....... (15) 


This equation is valid only for the case where the refractive 
index is the same as that for the medium in which it is immersed. 
that is, # = 1 fora droplet in a gas. 

Simpson?) proposed a solution for the absorptivity of a 
droplet when one takes into account refraction within the drop, 
but not external reflection. The resulting equation is: 


a, = m*/K?*r(2) ¢ 


X 


f.—2Kr [(m? —1) “) foe fm—1\" VT 
e 2K — £eS..6 
L L : m? Md J J 0) 


If one considers reflection and refraction of the radiant beam 
at every liquid-gas interface, (see Figure 1(b)) then a more 
complicated expression results: 


(2Kr + 1)e~*%r | + | 2K%r?/m? — 





pes 1 — e-2Kr cosy 


a, = 1 (1 _ Re)sin 20 ee Rp e728" cong dé ° ei) 


which must be solved by graphical methods. The parameters in 
this equation were previously defined by Equations (13) and 
(14). This expression resulted from an analy sis by geometrical 
optics or ray tracing, by Thomas‘! and Simpson‘), It should 
be pointed out that the absorptivity calculated trom Equation (17) 


| Reflection- refraction not 

a Considered- maQ 

|! Refraction only considered 
m=|.33 


ui Reflection-retroction 
Considered m«/33 





0 10 20 
Kr 3x0 


Figure 1(a) — Diagrammatic sketch of drop showing re- 
fracted and reflected radiation. 


Figure 1(b)—Plot shows effect of liquid absorption co- 
efficient and drop radius on droplet absorptivity. 
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is not entirely correct. In its derivation the phase shifts and lags 
were not considered and thus interference effects have not been 
included. However the absorptivity calculated from this equation 
is considered to be a good approximation to the actual absorptiv- 
ity, as long as the drop is fairly large. 

Absorptivities calculated from these equations for a real 
refractive index, # = 1.33 (this is the average real part of the 
refractive index of water for wave-lengths up to 15 microns) 
are shown in Figure 1(b). It can be seen that the absorptivity 
approaches an asymptotic value at large values of the parameter 
Kr. Thus the absorptivity of a drop may approach that of a 
black body (here a, ~ 0.95) if the drop radius and/or the 
absorption coefficient of the liquid is large.* As expected 
Figure 1 (b) shows that reflection and re fraction have a significant 
effect on the absorptivity of a drop. Certainly the addition of 
foreign solid material into the droplet will increase its absorptiv- 
ity above that given by Figure 1(b). 


Simpson“* also goes on to show that refraction greatly 
affects the distribution of radiant energy absorption within the 
droplet, for all radiant beams must pass through a relative radius 
r’ = r/m. Under certain conditions of absorption, this phenom- 
enon may cause a higher absorption within the sphere at this 


point than at the surface, and internal boiling may result. 


Returning to the cloud of droplets, consideration may be 
restricted for the present to the case where Kr is large enough, 
so that almost complete absorption does take place within the 
drop — at least to such an extent that reflection/refraction 
effects may be neglected. Then the exponent in Equation (7) 
would be of the form kJ, where kJ = a, ar,? N,. This is 
convenient, for then the methods described for gas radiation 
would be applicable. It is unlikely that the absorptivity of the 
spray droplets in this process will ‘be close to one, which would 
require relatively large values of K or 1, unless, of course, any 


solute in the drop greatly affects the absorption coefficient of 


the drops. Thus Mie scattering or reflection/refraction effects 
should be considered. Since, however, the complexities arising 


in considering these effects are enormous (and even more so if 
account is taken of the absorbing medium in an enclosure of 


finite dimensions), it appears advisable to use Equation (7), no 
matter what the absorptivity of the particles may be and to 
ascertain experimentally the error involved in making this 
simplification. ‘The deviation will arise as a result of the 
scattered radiation, which may be of greater importance when 
the particles are in an absorbing gaseous environment, At least 
the variables to be considered are identified; qualitative informa- 
tion is made available; the shortcomings of the simplification 
are known; and a starting point is established from this analysis. 


3. Cloud of Particles in an Absorbing Medium 


As suggested above, the attenuation of radiation during its 
passage through a cloud of particles can be expressed as a 
function of k. = Qyeye. wr,* N,. On the other hand, the attenuation 
by an absorbing gas is a function of the gas absorption coefhicient, 
k,. Each of these factors enters as an exponent for the exponential 
function. 


When both particles and absorbing gas are to be considered 
together in the same enclosure, these “resistances’’ operate in 
parallel and it can easily be dee that the exponents &, and ke 
should be added to give the attenuation of the traversing radia- 
tion. Thus over fraction ‘a’ of the spectrum the equation 
expressing this decrease in intensity is 


fig or (18) 


*Simpson’* determined a value of K 15.5 cm.~! for heavy oil fuel 
droplets. Thomas”* feels that a value greater than 100 cm.-' should 
be used for water droplets. Measurements similar to Simpson’s have not 
been made for water, although the absorption measurements by Brown, 
Sato and Sage” on thin water films show the absorption coefficient to 
be quite high. 
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and over (1 — a) of the spectrum, where only the cloud is the 
active absorber:* 


or for the full spectrum, the intensity after traversing a distance 
/ is given by: 
I=I,| ae 


—(& + ke)l + (1 = a)e~ ee (20) 


Some of the problems inherent in evaluating & will be 
discussed subsequently. 


4. Choice of Zone Size 


In order to employ the Hottel-Cohen method as outlined 
in section 1, it is necessary to choose zones small enough that 
they may be considered isothermal, and such that they have 
substantially the same view of all other zones. This latter 
restriction is necessary only for the proper distribution of 
reflective fluxes; in industrial enclosures, where surface absorp- 
tivities are relatively high, this factor is only of minor import- 
ance. It will be realized that, for this method to be rigorous 
(that is, the above stipulations rigorously adhered to), the 
zones must necessarily be differential in size. Differential 
zones are mathematically impossible, however, because of the 
integral character of radiation, and thus zones of finite size 
must be considered. Maximum accuracy is possible if the zone 
size is kept small; this increases the number of zones in any 
given cylindrical enclosure, however, and makes the calculations 
exceedingly long and laborious. 


A cylindrical zone of height equal to diameter may therefore 
be considered to be a satisfactory compromise. This implies 
immediately that radial temperature gradients are small, and 
that longitudinal temperature gradients within the enclosure are 
not excessive. The former restriction seems reasonable, 
least for columns of relatively small diameter. Also, a suitable 
average temperature for each zone may be chosen. 

Appendix illustrates the formulation of the equations 
governing the reception factors for this case. The derivation 
of these equations is similar to that presented by Hottel and 
Cohen“, It can be seen that these equations are based on the 
fact that a bundle of rays leaving a radiator proceeds in a straight 
line. As pointed out in the previous sections, some “diffusion” 
of radiation from this pencil of rays takes place because of 
particle scatter. In computing these view factors this is not 
accounted for; however an exponential law will describe the 
decay of radiation passing through a medium such as that under 
consideration — the problem is to ascertain the exponent. 


5. Energy and Material Balances on Each Zone 


The Hottel-Cohen method necessitates the formulation of 
heat and material balances on each zone in the system. ['xpressing 
the balance on G;: 








(Thermal radiation 
} into G,; from Gy, Ge | 
Jand S;, Ss. (direc " 
plus reflected ) 


from G; to all other 


(Thermal radiation 
) , 
(gas and surface zones, 


Enthalpy flux out of 
G; as sensible heat in 
| the flowing material J 


Enthalpy flux into G;) 
+< as sensible heat in the } 
flowing material J 


(Convective/conduc 
+. tive heat into G; from | 
adjace nt zones 


{Conduc tive /convec- 
+< tive heat flow out of > 
LG; to adjacent zones 


(Heat into G; as a 
+¢ result of chemical 
reaction 


f leat load within zone 
+ JG; due to chemical 
| reaction, evaporation 

etc, 








*This assumes that the droplets absorb equally over all wave lengths. 
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In the resulting series of equations, the wall temperatures 
are known but the longitudinal gas-temperature distribution is 
not. Furthermore, since the evaporative load in each zone 
depends on the temperature of the gas in that zone, this factor 
emerges as an unknown as well for each gas zone. The result 
is a very complex trial and error solution, the steps for which 
are briefly as follows: 
(1) Assume a temperature distribution; 

(2) From the drop-size distribution and drop velocity distri- 
bution at the top — or inlet — of the column and using (1), 
calculate the evaporation of a size interval (say, the interval 
from 30 to 40 microns, with d,ye = 35) with distance. Repeat 
for all size intervals. 

(3) From (2) ascertain the evaporative load in each zone; 

(4) Using (3) calculate the temperature in each zone from 
the heat and material balance equations; 

(5) Repeat steps (2), (3), ( 
bution. 

The details of these calculations will be the subject of 
subsequent papers. 


(4) for this temperature distri- 


On the basis of the balance equations, the mechanism of 
heat transfer in the column can be determined and this will 
help in estimating the temperature variation down the column. 

In that part of the column where evaporation takes place, a 
dynamic temperature balance in the gas is maintained. The 
heat transfer by radiation and convection from the walls and 
other gas zones to the zone in question must be the same as 
that which is transferred by convection-conduction from the 
gas to the droplets to effect ev aporation and to heat the generated 
vapours to the temperature of the gas in that zone. * The temper- 
ature of the gas will remain relativ ely constant** over a length 
of the column, as can be deduced from the following discussion. 
One must remember that the evaporation of a drop of diameter 
d; by convection/conduction effects is given by: 

(dm/d); = [hi Aa; (te — ta;))/X..... 2.20... (21) 


where ty is the temperature of the drop and will be very near 
the boiling point of the liquid. The apparent heat transfer 
coefficient, h;, for a drop of oe d; moving at a low relative 
droplet Reynolds number (16, 17) can be satisfactorily approxi- 
mated by the Ranz and Marshal Equation (15), corrected for 
the outward flow of cold vapor from the drop“'®; 


h; d; » In (1 + B) 


= 2 0.60 Re'!? 
i =~ (2+ 


PP sie sie ners (22) 
where B is Spalding’s Transfer number, ¢p(tg ta)/X (18). 
Thus an increase in gas temperature fg will decrease the amount 
of heat transferred from the walls to the gas. Since the resistance 
to heat transfer at the drop surface is low — that is /; is high 
relatively low temperatures will be observed in the upper zones 
of the tower. The governing step in this transfer process will 
be the heat flux from the wall to the gas. Since the gas throughput 
is small, it is expected that the main transfer process will be that 
of thermal radiation. The gas temperature in any given system 
and at any given wall temperature will be a function of the drop 
surface area available for convective heat transfer from the gas 
to the droplets. Thus, only when evaporation of the smaller 
of the size groups Is completed will the temperature of the 
increase. ‘This particular phenomenon will be discussed i 
greater detail in a subsequent paper. 


If the feed solution does not contain any solute or suspended 
solid then the drops will disappe ar upon evaporation. If a solid 
residue remains suspended in the gas, Nusselt“ and more 
tecently Sleicher and Churchill®® in a more rigorous mathe- 
matical development, have shown that these particles, subjected 
toa radiant flux, increase in temperature very rapidly and finally 
exceed the temperature of the gas by an amount dependent upon 


*Radiation from the gas to the droplets will be small in comparison to 
that transferred by the convection conduction mechanism. 

**Otherwise some of the heat transferred from the walls will be required 
to increase the sensible heat content of the gas. 
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their size, absorptivity and position in the cloud from the radiat- 
ing source. The temperature of the particle will be such that 
the amount of radiant heat it intercepts and absorbs, balances 
the amount of heat it transfers to the gas by conduction/convec- 
tion. Because of the magnitude of the latter effects, however, the 
temperature difference between particle and gas will be relativ ely 
small», The presence of particles thus should make the 
temperature of the gas increase at a faster rate than if they 
were absent. 


Conclusions 

In the analysis which has just been presented, an attempt 
has been made to develop a method of solution capable of 
yielding results sufficiently accurate for engineering design 
purposes. Although most of the assumptions which have been 
made appear reasonable, it must be recognized that the proposed 
approach is incapable of dealing with ‘the extremely complex 
effects resulting from the scattering, reflection and refraction of 
radiation by the particles in an ‘absorbing medium of finite 
dimensions. The use of Equation (7) which is recommended 
in this study may prove to be an over- simplification, but appears 
to be a logical starting point for the present.* 

Before the method of solution can be used, it is necessary to 
know the reception factors for a cylindrical system containing 
an absorbing medium. These have been calculated (see Appendix 
A) and the results are shown in Figures 2 to 8. 

To estimate the temperature profile i in the column by means 
of the heat balances discussed in section 5, it is necessary to 
calculate the magnitude and distribution of the radiant heat flux. 
This in turn necessitates the calculation of total interchange 
areas between zones, making due allowance for the multiple 
reflections in the system. These total interchange factors can 
be built up from the reception factors, as shown in Appendix B. 
It is interesting to note that the latter procedure also permits 
the estimation of the mean beam length for any cylindrical 
system, as discussed in Appendix C. 


— -_ — > 


r- 0.01 


- 0.005 





ee 0.0001 
0 Os 


10 1S 20 25 30 


ks 


Figure 2— Reception factors for a cylindrical column 
containing an absorbing medium /,,,, Vicia Dean 


*In a recent analysis, presented at the New York meeting of the A.S.M.E., 
Dec. 1960, Hottel and Sarofim have indicated how the zoning technique 
described in this paper could be modified to account, in an appropriate 
manner, for the scattering of radiation in particulate systems. 
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Figure 3— Reception factors for a cylindrical column 
containing an absorbing medium f,,., f,.,5 fs. 


One remaining difficulty should be mentioned. Because 
evaporation of the droplets in the column necessarily means a 
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Figure 5 — Reception factors for a cylindrical column 
containing an absorbing medium @ = 1 — f,,.. 
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Figure 4— Reception factors for a cylindrical column 
containing an absorbing medium f,,,, f..,5 Fees 


reduction in their size of their complete disappearance, the 
absorption coefficient &. for the cloud is constantly decreasing 
with distance down the column. Approximate methods are 
available to deal with this problem which are discussed in 
Appendix D. 

The details of an experimental study in which the above 
method of solution has been successfully applied will be presented 
in a subsequent paper. 
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Figure 6 — Reception factors for a cylindrical column 
containing an absorbing medium f,,,, f,.,5 foc. 
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Figure 7— Reception factors for a cylindrical column 
containing an absorbing medium f,,,, f,,,5 f4s5° 
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Figure 9—Differential elements upon which the formula- 
tion of the reception factors is based. 
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Nomenclature 

Roman Symbols 

A, = surface area, ft.? 

ay = weighting factors, dimensionless, (a in Equations 


(3) and C3) 
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Figure 10(a)—Build-up of area interchange factors. 


Figure 10(b)—Mean beam length for an n-zone system. 


An, On = complex functions of x, defined in Reference (4), 
dimensionless 
b = side of cube or square 
B = Transfer Number: c,AT/A, dimensionless 
Cy = heat capacity, (B.T.U.)/(Ib. (°F.) 
d = diameter of droplet, ft. 
Ei(—y) = exponential integral of y 
f = reception factor between unit zones in a black- 


walled cylindrical system containing an absorbing 
medium — denotes the fraction of radiation leaving 
one zone (first letter) which reaches and is absorbed 
by another (second letter) 


(1) f,, = between two vertical surface zones 

(2) f., = between an end surface and a vertical 
surtace zone. 

(3) fee = between two end surface zones 

(4) fy, = between a vertical surface and a gas zone 

(5) feg = between an end surface and a gas zone 

(6) f,, = between two gas zones 


(subscripts defined on Figures 2 to 9), dimensionless 
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Fiz = direct geometric view factor between two surfaces 
in a black-walled system containing a non-absorbing 
medium (see Reference 2), dimensionless 

FF = total interchange factor in a grey-walled system 
containing a non-absorbing medium (see Reference 
2), dimensionless 


F,,, Fe = — elemental reception factors corresponding to 
ffactors, but between differential elements dh high 

g,G = gas zone with numerical subscript to designate the 
particular zone under consideration 

ze = direct interchange area between any two gas zones 
in a black-walled system, ft.? 

gs = direct interchange area between a gas and surface 
zones in a black-walled system, ft.? 

GG = total interchange area for a grey-walled system 


containing an absorbing medium which regulates 
the amount of radiation which reaches any one gas 
zone due to emission by another, including reflection 
at all grey surface zones, ft.? 

GS = SG = total interchange area for a grey-walled system 

containing an absorbing medium which regulates 

the total interchange between a gas and a surface 
zone, including reflections from other surface zones, 
<* 

vertical separating distance, ft. 

heat transfer coefficient for drop of diameter d; 

intensity of radiation from a surface 

intensity of radiation from a surface at ] = 0 

subscripts to represent any gas or surface zone 

thermal conductivity of gas, (B.t.u.)/(hr. (ft. )(°F.) 
absorption coefficient, ft. 

absorption coefficient of cloud, ft 

absorption coefficient of gas, ft.- 

absorption coefficient of liquid, ‘ 

thickness of absorbing medium, ft. 

mean beam-length, ft. 

refractive index made up of a real and an imaginary 

part 

integer in E quations (9) and (10) 

real part of refractive index 

r imaginary part of refractive index 

number of particles in size group 1 

any whole number, including zero 

number of particles per unit volume with radii in 

the interval dr 

N, = particle concentration, (number of particles)/(cu. 

ft.) 

= heat transfer rate, (B.t.u.)/(hr.) 

efficiency factor for extinction in the forward direc- 

tion, dimensionless 

radius of droplet, ft. 

radius of droplet in size group 1, ft. 

surface mean radius, ft. (defined in Equation (12)) 

fraction reflected at ‘surface, Equation (13) 

real part of expression in Equation (10) 

subscript, refers to radiant heat 

scattering cross-sectional area, ft.? 

extinction cross-sectional area (absorbed plus scat- 

tered), ft.? 

Si, Se = refers to surface zone, subscript number designates 
specific one 

35 = direct interchange area between any two surface 
zones in a black-walled system containing an 
absorbing medium, ft.? 

ve = direct interchange area between a vertical surface 
and an end surface in a black-walled system contain- 
ing an absorbing medium, ft.? 

3Z = direct interchange area between a surface and a 
gas zone in a black-walled system containing an 

ns absorbing medium, ft.? 

SS, SG = total interchange areas for a grey-walled system, 
which regulates the amount of radiation which 
reaches a surface or gas zone, respectively, due to 
emission by a surface zone, including reflections at 
all grey surface zones, ft.* 

& = absolute temperature, °R. 
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te = gas temperature, °C. 
ta = particle temperature, °C. 
V = gas volume, ft.é 
7 a ‘ ac , > f 3 
t = total gas volume, ft. 
x = parameter: 2mr/X\, dimensionless 
Z = characteristic length of an enclosure, ft. 
7, 0,2 = cylindrical co-ordinates 


Greek Symbols 


a, = particle absorptivity 

Y = extinction coefficient, ft.~} 
6 = cylinder diameter, ft. 

€c = gas emissivity 
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€x, 22 = gas emissivity evaluated at one mean beam-length 
and at two mean beam-lengths, respectively 


6 = time, hours 
6 = incident angle, in Equations (13) and (14) 
A = wavelength, ft. 
A = latent heat of vaporization, in Equation (21) 
y = refracted angle, in Equations (13) and (14) 
V = factor by which Z must be multiplied to yield the 
mean beam-length 
® = escape factor = 1 — fer, 
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APPENDIX A* 


Derivation of Reception Factors for a Cylindrical 
System 


This discussion is restricted to a cylindrical enclosure contain- 
ing an absorbing gas of constant absorption coefficient. This 
enclosure is divided into a number of surface and volume zones, 
each of height equal to the diameter of the cylinder (hereafter 
referred to as unit-zones). The problem is to find the reception 
factors for this system, defined as the fraction of the total 
radiation that leaves any radiating zone (surface or gas) which 
is intercepted by another zone, taking into account the absorption 
of radiation by the intervening medium. ‘Thus this factor 
corresponds to the direct-view factor, F, previously used by 
Hottel », but now for a system containing an absorbing medium. 


*The complete derivations contained in Appendices A and D have been 
deposited as Document No. 6842 with the ADI Auxiliary Publications 
Project, Photoduplication Service, Library of Congress, Washington 25, 
D.C. A copy may be obtained by citing the Document No. and by 
remitting $3.75 for photoprints, or $2.00 for 35 mm. microfilm. Advance 
payment is required. Make cheques or money orders payable to: Chief, 
Photoduplication Service, Library of Congress. 
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1. Surface-to-Surface Interchange. 


By definition, the one-way radiation exchange between any 
two black- walled. surface zones is given by:** 


Is} = A,, ° Tess . oT,,* 


In terms of differential areas (Figure 9), this one-way exchange 
is also known to be directly proportional to the area of each as 
viewed from the other (Lambert's Cosine Law), and inversely 
proportional to the square of the separating distance R. In this 
case the intensity at dA» will be further reduced, for account 
must be taken of the transmissivity of the gas over the distance R. 
Expressing this mathematically: 


dijy-s9 = (1/t) (dA cos 8; dA» cos /R2) (e #®) (6 T,,4) . .(A2) 


When this expression is integrated over the two finite areas, 
Equations (A1) and (A2) are equal and thus f,,,. can be evaluated. 
For the general case of zone i and 7 = i + N (where N is any 
whole number including zero) and when all length terms are 
normalized with respect to the cylinder diameter 6: 


or 1 N+1 
fisyy = a dd, | dd» d(22/6) d(z,/8) 
‘ 16m? P . 


0 o N 


— cos(¢: — dri)]*. « 
where 


(R/6b)? = sin? (do — gi)/2 + [(22 


he nomenclature adopted to describe f-factors is shown on 
each f-plot, (see Figures 2 to 8). The final letter subsc ript 
indicates the receiving unit-zone, and its type, just as the first 
letter shows the primary zone. The number on the final subscript 
indicates the number of unit zones that the receiving zone is 
removed from the primary one. In each case N is this number. 


2. Gas-to-Gas Interchange. 


Similarly, the generalized view factor for gas zone i to gas 
zone 7 = 7 + N becomes: 


kd 1/2 
ley = AK « dd» [ow 6)d(r,/8) (72/6 )d(r2/8) 
= o 
N+1 


d(z1/5) —k6(R/6) 


° N 


d( 22/5) e /(R/6)? ....(A14) 


3. Surface-to-Gas Interchange. 


For the paws case of surface 7 to gas volume 7 = i + N 
or 1 N+1 
: = ° ° JN 
‘es »/5)d 6) [1 (t= (12 6)cos(ds = oi ] 


Direct analytical integration of these integrals was not 
possible. They could, however, be partially solved analyticaily 
(the complete mathematical analysis is presented in Reference 
24). The final integrations were evaluated numerically, and the 
f-factor working plots are shown, together with the nomenclature 
adopted, on Figures 2 to 8. 

The numerical integrations were not carried out for zones 
greater than two removed, since the view-factors for these zones 
were very small. Approxim ite techniques can therefore be 
employed without causing serious error to the overall result. If 
the separating distance between unit-zones becomes large with 
respect to the diameter of the column, then an average trans- 
missivity can be used for all points in either zone. It can be 
shown that the equations for F,,, F,, and F,, at 86 = 0 are special 
cases of Equations (A25), (A26) and (A27) respec tively. 


1 h(2h? + 3) 
F,, (h, R86 = 0) = =| 2 — — ; (A31) 
2 (h? + 1)3/2 
2h? + 1 
F., (h, k6 = 0) = 2 oe - 2h (A32) 
(1 $ h?)!' 2 


*°Simpson“ has shown that this equation is only valid when the re- 
fractive index is one or very close to one. If it is m, then @T)4 should 
be multiplied by m*. This correction will be omitted in this discussion, 
only because radiation in a gaseous medium is being considered here; 
for other materials the refractive index must be considered. 
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Fre (h, RO = = 2 [h? — h(1 + h?)”? + 1/2]... .(A33) 


These equations have been integrated with respect to h to give 
fug = 0.004554, fug = 0.01114, fg = 0.01515. To find the 


values of these factors at any &6, it is recommended that these 


( 
values by multiplied by e BH R./d) where R, is the centre-to- 
centre distance. Since this procedure is not possible for the 
surface-gas and gas-gas factors, it is recommended that they be 
found by proper addition and subtraction of the surface-surface 
factors at the required ké. For k6 > 1, these factors are insignifi- 
cant and can be omitted. 


APPENDIX B 


Calculation of Direct-Interchange Areas 


So far, evlindrical surface and gas zones of a size equal to a 
right circular cylinder have been considered. The one-way 
radiation is ~ n given by equations such as (A1). 

In some cases, however, the conditions may be such that 
larger lecchermal zones may be postulated. In these cases, it is 
possible to build up interchange factors from the existing plots. 
Using Hottel and Cohen’s nomenclature, as discussed in their 
paper"), one may define one-way radiant flux for surface-to- 
surface, surface-to-gas, gas-to-surface, and gas-to-gas respectively 


as: 


Ysjsy = (Si5;) rl PR pete oy .(B1) 
Qsiy = (81%) OT,,*. . . .(B2) 
Ge; = (2i5;) OT,,*. . .. (B3) 
ee; = (21%) OT,,*...« (B4) 





Th» factors s;:,, ;2 si, gig; have the dimensions of area and are 
defined as direct-interc hange areas the first letter with its 
subscript indicating the sender and the second the receiver. 
These factors apply only in a black-walled system and are 
analogous to the factors AF previously used in a black system 
containing a non-absorbing medium. Just as in that case: 


Dis s; + rs = A,, : ak (BS) 


J 


i 


j 


and 


a eee 


It also follows from the application of Kirchhoff’s Law to an 
isothermal system that 


OS. = Sift... 6.5 : o<e Ge 
or by definition: 42 V6, . fers, = Asy - Sere, + - 5 ccteta een (B8) 
or for a right circular cylindrical zone: 26f,,5, = faye. - - - - --(B9) 
Similarly gico = geg; and s\s2 = S2 ‘Ss, Which can also be seen from 


the symmetry of the defining equations for the f-factors. 


The remaining problem is to find a factor such as gig, for a 
cylindrical-zone larger than the basic one used to formulate the 





f-fac tors. To find the direct- interchange area gig 1 for a gas zone 
comprising two right circular cvlinders as shown in Figure 10(a), 
the following formulation is necessary: 


“18 — 153 A i . . es 
z 1 = 9% = 4k [Peete + Sassy + Sects + Seana] ----- -( B10) 
ao ‘1 


where the first two f-factors are of the fy, type, while the last 
two are of the f,,, type. 


gsy 


\ similar expression results for sys). 


To find $\ Se 


sj 59 cnsinias 


6S = wo? (f. « Be Pes) 


(Bil) 


where the fust factor is of the fs, type, while the second is of the 


fss, type. 


APPENDIX C 


Estimation of Mean Beam Length 


By means of the foregoing procedure it is possible to determine 
the mean beam length for any cylindrical system. 
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The mean beam length of an enclosure containing an absorbing 
or emitting medium is defined as the radius of a hemisphere 
containing the same medium at the same temperature which will 
cause the same radiant flux density (B.t.u./ft.? hr.) on a differ- 
ential area at its centre as that which the walls of the enclosure 
experience from the total gas it contains. Thus, what is required 
is the flux density leaving the gas mass. Referring to Figure 10(b), 
the total cylindrical volume V; is made up of n-basic right-circular 
cylindrical unit-zones. The flux leaving this gas mass is 


4k V,®,, = 4kV; — S6in<ss CT) 
It can easily be shown that 


Zig = 4(RO) (7162/4) 


[fguq + 2(m — 1) fee, + 2(m — 2)fery + 2(m — 3)feegt+...] . (C2) 


The flux at the walls resulting from the contained gas is 
a(4kV®,,) oTy,!.... (C3) 


and if this is divided by the total area of the cylinder (including 
the ends), the radiant flux density can be calculated. The 
weighting factor ‘a’ necessary for this equation is valid only over 
part ‘a’ of the spectrum. The emission of the gas in a hemisphere 
per unit of surface at the centre of a hemisphere of gas of radius 


L is: 


éc, oTy, = a(t —e"*)oTys.... (C4) 
Equating (C3) and (C4) 
_pys  4kV Oy 
l—e kVO _ ——.. ACS) 
Aig 


where WV is the factor by which the characteristic dimension of 
the system must be multiplied to give the mean beam length of 
the system. Since the evaluation of ®y, involves a knowledge of 
&6 and thus of L, Equation (CS) must be solved by trial and error. 
WV versus £6 for a cylinder (length = 106 and 156) is shown in 
Table 3. Only a small variation in V is observed. The value 
recommended by Port (23) from his studies on rectangular 
parallelepipeds is included in this Table —i.e., L = 0.85 (4V/A) 
for average values of Pcl. For infinite cylinders 0.906 is rec- 
ommended®). 


APPENDIX D* 
Treatment of Varying Absorption Strength 


Since the temperature of the gas varies and its emission or 
absorption is known to depend on temperature, it follows that 
the absorption coefficient, as well as the weighting factor ‘a’, 


*See footnote for Appendix A. 
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will vary as well. Moreover, a more serious variation in absorption 
coefficient of the cloud occurs because 

(1) the size of the drops decreases as evaporation proceeds; 

(2) as evaporation proceeds more vapor is evolved and thus 
the voidage increases or the number of particles per unit volume 
decreases; 

(3) the voidage increases as a result of the increased temper- 
ature, increasing the volume of unit mass of gas; 

(4) in some cases the droplets may disappear altogether; 

(5) the particles may explode — this explosion may occur 
when treating some solutions where an impervious solid shell 
forms on the surface of the drop and the vapor generated 
internally causes an increased pressure inside the droplet, which 
finally overcomes the forces holding it in. 

Theoretically it is possible to account for this changing 
absorption strength down the column by starting with the 
numerical values of F, and in the final numerical integration 
with respect to (2/6), (22/6) and/or (22 — 2)/6, by using an 
integrated average &6 over the separating distances in question, 
Vi1Z. 


(R8 )ave = b6 dz ie es. sven hanwe (D1) 


However the labour involved would. be prohibitive and thus 
again a simplification is necessary. It can be shown* that the 
vaiious direct-interchange areas and f-factors can be corrected 
for variable absorption strength approximately in the following 
manner: 








CRY Ge GR in I oink ci cid casero (D7) 

) == . (26) gas zone -— 

(2) (sg or fre )var 25 = (hee (sg or fg) ave abs (D8) 
_ kb) zone | 

(3) (gs or fysjvar 25 = ( i. (gs or fis.aye ae (D9, 
aa (B6)c, (2 )eg — 

(4) (g)var es = 


—— ie ara (D10) 


Sample calculations illustrating the use of these equations will be 
presented in subsequent papers. 


°See footnote for appendix A. 
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An Apparatus for the Measurements of 
Light Transmittance of Emulsions’ 


J. A. CENGEL’*, J. G. KNUDSEN’, A. LANDSBERG: 
and A. A. FARUQUE 


This paper describes an apparatus for the meas- 
urement of the light transmittance of emulsions. A 
photoelectric technique was employed to determine 
the effect of flow rate, mixing time, and dispersed 
phase concentration on the light transmittance of a 
liquid-liquid emulsion composed of a commercial 
petroleum solvent and water. For a given concentra- 
tion, light transmittance was independent of flow rate 
and became constant after 10 minutes of mixing in 
the test apparatus. Light transmittance decreased with 
solvent (dispersed phase) concentration. Up to 5% 
solvent concentration, the light transmittance-concen- 
tration curve could be used to measure solvent con- 
centration accurately. 


he flow characteristics of liquid-liquid emulsions are of inter- 

est in many processes in which two liquid phases contact each 
other. The rheological properties of such emulsions depend 
upon such factors as the proportion of the two liquids present, 
the interfacial area, the type and amount of emulsion stabilizer 
present, and the size of the dispersed phase particles. In addition, 
experimental values of the rheological properties are dependent 
upon conditions under which measurements are made. 

This paper describes a photoelectric technique for measure- 
ment of the light transmittance of liquid-liquid emulsions and 
demonstrates the effect of dispersed phase concentration, emul- 
sion flow rate, and mixing time. The emulsion studied was 
made up of a commercial petroleum solvent as the dispersed 
phase with water as the continuous phase. This study is part 
of an over-all research program on the flow and convection heat 
transfer characteristics of unstable liquid-liquid dispersions. 
The various properties (including viscosity, the thermal con- 
ductivity, light transmittance, and interfacial area) are deter- 
mined and related to dispersed phase concentration and to the 
momentum and heat transfer characteristics of the phase fluid. 

Langloise and Gullberg™ used a light transmittance technique 
to determine the interfacial area in a two-phase liquid-liquid 
emulsion. These workers, however, studied emulsions in which 
the dispersed phase droplets were relatively large and developed 
a relationship by which light transmittance was related to inter- 
facial area per unit volume as follows: 


Felt BA WP Vesacwexiaers (1) 


SRR SR SESE SSR ESSER TEESE SEE SEES EEEEREE SEAHEC REERES ERE EERE REE E REE REE ee Eee ee, 


‘Manuscript received March 17; accepted June 28, 1961. 

*Department of Chemical Engineering, Purdue University, Lafayette, 
Indiana, U.S.A. 

‘Department of Chemical Engineering, Oregon State College, Corvallis, 
Oregon, U.S.A. 

Contribution from the Department of Chemical Engineering, Oregon State 
College, Corvallis, Oregon, U.S.A. 
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They considered that the constant B was independent of 
volume fraction of dispersed phase. In_the present work, 
however, it was not possible to determine the interfacial area 
of the emulsions studied since the dispersed phase particles 
were very small and light scattering brought about by these 
particles prevents Equation (1) from being valid. ‘As the 
dispersed phase particles become smaller and smaller, the light 
transmittance becomes independent of interfacial area for a 
given concentration. At small dispersed phase concentrations 
light transmittance is a sensitive function of concentration, 
howev er, and may be used to predict concentration. 


Equipment and Procedure 

The experimental system is shown in Figure 1. The photo- 
electric apparatus consisted simply of a constant voltage light 
source and a photocell. Figure 2 is a detailed drawing of the 
light and photocell probes. 

The light-source probe (8) was mounted on the main piping 
system (16), perpendicular to the flow, by soldering a brass 
fitting (14) into a 5¢-in. hole in the pipe wall. The piping 
system and the light-source probe were sealed from one another 
by the glass window (15) in the — steel light directing 
tube (9) and a packing gland (13), which was forced into the 
stuffing box by the fitting (12). The light-source probe (8) 
was soldered to a ring (10) and this combination was held to 
(12) by three brass screws (11). The end of the light-source 
probe was closed by a micarta end-piece (3) held in place by 
binding post (2), which also served as a ground connection. 
Two binding posts (1) and three adjustment screws (4) were 
fitted into the end piece. The aluminum lamp base (6) and the 
lucite holder (5) are adjustable to give the proper illumination 
from the lamp. 

The photocell probe was soldered to the main piping system 
directly opposite the light source tube by means of fitting (17) 
in a 114-in. hole. This probe was sealed from the system by 
window (20) in the photocell supporting tube (24). The packing 
was held in place by gland (18) which was forced into the stuffing 
box by fitting (19). The photocell fitted into the socket mounted 
in lucite (21) and was attached to the micarta end piece (22). 
Binding posts (26) supplyi ing power to the photocell were 
mounted on the aa piece. The entire photocell mounting was 
held in place by set screw (25). 

A Delco 6-volt lead storage battery was the voltage source 
of the light bulb (7). The current was first directed into an 
exterior electrical system so that a specified voltage (usually 4.5 
volts) could be maintained at the light bulb. ‘The voltage 
delivered across the light bulb was checked before each reading 
to insure that all data were taken under identical conditions. 

Che transmitted light received by the photocell tube (23) 
was converted into a potential which was measured by a null- 
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Figure 1—Experimental equipment. 


point potentiometer. A Leeds and Northrup Model 2430 
galvanometer used to observe deflection was external to the 
potentiometer system. The face of the galvanometer was 
calibrated from 0 to 100 in increments of one, to allow the 
estimation of percent changes in light transmittance. The 
instrument was set to read zero with the light-source off and 
100 with the light source on and with clear water flowing in 
the main piping system. Sensitivity of the galvanometer was 
+1%. 

The electrical system is schematically shown in Figure 3. 
The symbols represented are as follows: 


Bl 90-volt battery (RCA VS090) 

B2  6-volt lead storage battery (Delco dry charge) 

B3 4 mercury cells (Mallory ZM-9) 

C1 Two sets of contacts for phototube (RCA 1P41) 

2. Two sets of contacts for light (GE 82, 6-volt) 

C3 Galvanometer connections (Leeds & Northrup 2430a) 

Ri Coarse adjustment rheostat (10 turn 20,000 ohm Helipot) 

R2 Fine adjustment rheostat (10 turn 20,000 ohm Heliport) 

R3 Load resistor | megohm) 

R4 Coarse adjustment rheostat (5 ohm rheostat) 

RS Fine adjustment rheostat (10 turn 25 ohm Heliport) 

R6 Load resistor (50,000 ohms) 

R7 Balancing voltage set potentiometer (10 turn 50,000 ohm 
oo 

R8 Voltage resistor (10 ohm) 

R9 Sensitivity lowering resistor (50,000 ohms) 

R10 Sensitivity lowering resistor (1000 ohms) 

R11 Sensitivity lowering resistor (50 ohms) 

SI Double pole double throw circuit selector switch) 

$2 Single pole double throw push button 

$3 Single pole double throw cell selector switch 

$4 Double pole single throw push button 

S5 Single pole single throw light switch 

S6 5 position sensitive selector and galvanometer switch 


= 


The emulsion, composed of a commercial petroleum 
solvent, (Shellsoly 360) and water was prepared in a stainless 
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Figure 3—Wiring diagram. 
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Figure 2—Light and photocell probes. 


steel, jacketed tank using a variable speed stirrer for agitation, 
The solvent had a viscosity nearly equal to that of water in 
the temperature range studied. The emulsion was circulated 
through the piping system shown schematically in Figure 1. 
The system was raade of 114-in. standard brass pipe while the 
section containing the probes was 2-in. standard brass pipe. 

The supply tank and main piping system were flushed with 
solvent several times before runs were made. A calculated 
weight of solvent was added to previously weighed amount of 
water in the supply tank. To obtain thorough blending of the 
two liquids, mixing time was usually 2-3 hours. Agitation was 
continued during the test. 

The mixture had a milk-white appearance characteristic of 
many liquid-liquid emulsions. If the agitation was stopped, a 
clear layer of solvent immediately appeared at the surface of the 
liquid in the supply tank, indicating instability of the emulsion. 

Samples of the dispersion were analyzed periodically by 
allowing the emulsion to separate in a graduated cylinder. 


Experimental Data 

The light transmittance of mixtures containing 0.49, 0.97, 
1.92, 2.87, 3.82, 4.75, 19.4, 34.3, and 49.2 volume percent 
solvent was determined at probe separations of 4%-in. and 14-in. 
In addition, the effect of mixing time and mass flow rate on the 
light transmittance of the latter four mixtures was investigated. 
Details of the experimental program are available in Reference 
(2). 


Results 


Figure 4 is a plot of percent light transmitted versus length 
of mixing time. At zero time, the required amount of solvent 
is added to the mixing tank and circulation and agitation initiated. 
The percent transmittance decreases rapidly and becomes con- 
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Figure 4— Amount of light transmitted as a function of 
mixing time. 
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Figure 5—Effect of flow rate on amount of light transmitted. 
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Figure 7—Relation between solvent concentration and light 
transmittance. 


stant after about 10 minutes. It was observed by visual inspection 
and by measurement of friction losses that the emulsion was not 
yet thoroughly mixed and the dispersed phase particles were still 
being reduced in size. These results indicate that after about 
10 minutes of mixing, the dispersed phase droplets were sufti- 
ciently small that the transmittance was independent of droplet 
size and was unchanged by any further size reduction. 

Figure § is a plot of percent transmittance versus flow rate. 
There appears to be no effect of flow rate on percent trans- 
mittance indicating that any change in drop size or drop-size 
distribution could not be detected by the present apparatus. 


The effect of solvent concentration is shown in Figures 6 
and 7 in which percent transmittance is plotted versus volume 
percent of solvent in the emulsion. Figure 6 shows a rapid drop 
in percent transmittance with concentration in the range up to 
§% solvent. Beyond this the transmittance decreases more 
slowly with concentration ultimately becoming constant in the 
range 35 to 50% solvent. 
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Figure 6—Relation between solvent concentration and light 
transmittance. 


These results indicate, however, that the light transmittance 
technique would be useful in measuring low solvent concentra- 
tions once the light transmittance-concentration curve is estab- 
lished. Figure 7 shows the results obtained for two probe 
spacings up to 5% solvent. In this range the light transmitted 
is very sensitive to the concentration ‘changes. Beyond 5% 
solvent, however, the light transmittance becomes oulliciunite 
insensitive so that accurate measurements of concentration could 


not be made from light transmittance determinations. 

In summary, an apparatus is described to measure the light 
transmittance of unstable liquid-liquid — flowing in 
circular tubes. Experiments showed that for a given solvent 
(dispersed phase) concentration in water, light transmittance 
became constant after 10 minutes of mixing and was independent 
of mass flow rate in the pipe. E ight transmittance Is very sensi- 
tive to concentrations in the range up to 5% solvent. In this 
concentration range and for a given probe separation, the light 
transmittance-concentration curve may be used to measure 
solvent concentrations accurately. Beyond 5% solvent lght 
transmittance becomes insensitive to concentration. 

Work is presently under way to determine the droplet size 
of the various dispersions. At the low dispersed phase concen- 
trations where the light transmittance is a sensitive function of 
concentration studies are being made on the effect of flow rate 
and mixing time on surface area to determine the relation if any 
of these factors on the light transmittance of these dilute dis- 


pe r sions. 
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Nomenclature 

A Interfacial area per unit volume 

B Specifying constant dependent on the ratio of refractive 
indices 

J, = Light intensity incident to the emulsion 

I = Light intensity emergent from the emulsion 
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Optimum Performance of Combined Flow 
Reactors under Adiabatic Conditions’ 


A, CHOLETTE? and J. BLANCHET? 


Irreversible first- and second-order reactions are 
studied for combinations of tubular and CST reactors 
when endothermal or exothermal adiabatic conditions 
prevail. 

Typical results presented for MT and TM models 
show the variation of either the conversion or the 
relative residence time as a function of the mixing 
level of the system. 


A graphical method, applicable generally to MT 
reactors, provides a simplification over that of suc- 
cessive approximations required for model TM. 


For the conditions investigated broad limits are 
indicated of the optimum reactor combinations. 


gc studies in chemical reaction engineering indicate a 
growing trend toward optimization techniques and empha- 
size the importance of evaluating conditions leading to optimum 
performance®:?:*-*, 

Since flow reactor operation is highly susceptible to two 
important factors, temperature and the extent of mixing, each 
of these has been the subject of numerous investigations. Some 
have dealt, in particular, with temperature profiles in tubular 
reactors®:®.7-) or with temperature sequences in continuous- 
stirred-tank-reactors (CSTR). Others have been suggested 
by the extent of backmixing or partial mixing and resulted in 
such studies as on the distribution of residence times and axial 
mixing ™!.12,13,16) | 

In the case of flow reactors which do not exhibit true plug 
flow nor CSTR behavior the use of models is being resorted to 
more and more as a fruitful method of approach. Thus the 
performance of partially mixed reactors may be analysed through 
models consisting of combinations of tubular and well-mixed 
sections 415,16.17,18,19 | In some instances actual results have 
shown good agreement with the proposed models”. 

The behavior of some model reactors, for isothermal opera- 
tion, has been —, recently®!©, It involved the introduc- 
tion of a parameter e defining the level of mixing as that 
fraction of the total edie of the reactor which was considered 
to be perfectly mixed. Interesting results were thus obtained 
and more were anticipated when considering such models in 
non-isothermal operation. This has triggered ‘the present study 
on irreversible first- and second-order reactions carried out 
adiabatically in combinations of reactors. 


'Manuscript received December 7, 1960; accepted July 20, 1961. 
*Chemical Engineering Department, Laval University, Quebec, Que. 
Based on a paper presented to the C.1.C. Chemical Engineering Confer- 


ence, Quebec, Que., November 7-9, 1960 
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Les auteurs étudient les réactions irréversibles du 
premier et du deuxiéme ordre pour des combinaisons 
de réacteurs tubulaires et bien agités, dans le cas de 
réactions endothermiques ou exothermiques. 


Les résultats caractéristiques pour les modéles MT 
et TM montrent la variation de la conversion ou du 
temps de résidence relatif en fonction du niveau 
d’agitation du systéme. 

Une méthode graphique, qui s’applique de facon 
générale aux réacteurs MT, est plus simple que celle 
par approximations successives que l’on utilise pour 
les réacteurs TM. 


On indique, dans les grandes lignes, les combinai- 
sons de réacteurs qui donnent un rendement maxi- 
mum. 


ADIABATIC REACTORS 

Since many factors affect the rate of reaction, those con- 
sidered in this work have been limited to temperature and 
concentration, in ideal homogeneous systems. The analytical 
treatment of combinations of reactors requires that one consider 
first the overall effects of variations in temperature and concen- 
tration on the behavior of both tubular reactors and CSTR’S. 
This can be facilitated by expressing temperature effects in terms 
of variations in conversion. The latter, in first- and second-order 
reactions of the type: 

A — Products 
2A — “ 


A+B-— ‘ 


is, by convention, defined as the ratio of the quantity of A 
reacted over the original amount present. 

In any adiabatic system the heat effect can be related to the 
conversion by the following expression provided the heat of 
reaction AH, and the specific heat of the reaction mixture 
remain constant: 

T = T, + dx. . . 3 ove 
where 6, the maximum temperature variation for complete 
conversion, is a function of AH/z and Lijec ,;. 

For convenience Equation (1) may be rewritten in a dimen- 
sionless form: 


Similarly, the usual Arrhenius equation can be modified from 
hk = 2¢~2/RI - 
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TABLE 2 
DESIGN EQUATIONS FOR ADIABATIC FLOW REACTORS 


Order of . : a 
earth Tubular CSTR 





1. Xs dx x 
6, -| a Cig en 











o Al — x) Rk(1 — xr) 
2 CA pt = Gi Ou 
Xt dx Xu 
o Mi —x)(Rrp—x) | Mi—xw)(Re—xw) 
E 1 bx 
Rie r 
to: wee tC = he era wee eea (4) 


The latter can be further transformed after introducing 


Equation (2): 
E cee 
RT, % 


1 
k=ke a 


k, being the initial reaction rate constant, at T,. 

It can be seen from this expression that, all other conditions 
being equal, the reaction constant may be assumed to be a sole 
function of x. It then becomes easy to compare different kinds 
of reactions if the parameters EF, b, K and T,, or combinations 
of them, are properly chosen. 

Although, in practice, values of the above parameters may 
vary over a wide range, they have been limited in the present 
study to those specified in Table 1. This has been done for 
convenience on account of the long trial and error calculations 
involved due to the lack of proper computing facilities. Had 
these been available, an investigation of much wider scope would 
have been possible. 
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Figure 2—Variation of the conversion with residence time 
and initial concentration for second-order exothermal 
reactions. 
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Figure 1—Variation of the conversion with residence time 
for first-order exothermal reactions. These “working curves” 
are for conditions of reaction specified in Table 1. 


Some considerations may be in order at this point on tubular 
and CST reactors taken separately. Equations"'*.®) developed 
for these particular cases are presented in Table 2 for first- and 
second-order reactions, in constant density systems. 


Solution of these equations, —— for the variation of k 
and using the specifications of Table has resulted in 18 
working graphs, three of which are shown <A Figures 1, 2 and 3. 
The numerous curves are for first-order reactions and for second- 
order reactions with values of R; equal to 1.0, 1.25, 2, 5 and 10 
respectively. They were prepared, in each case, for isothermal, 
endothermal and exothermal conditions. 


It can be observed in Figures 1 and 2 that the exothermal 
case exhibits curves which behave differently from those for 
isothermal or endothermal conditions. At small and moderate 
conversions the CSTR curve lies over that for the tubular 
reactor but at higher conversions it reverts to the relative 
position exhibited in isothermal“® and endothermal reactions 
as in Figure 3. 

While a CSTR is at a disadvantage with respect to a tubular 
reactor from the point of view of a lower effective concentration, 
it takes advantage on the other hand in exothermal reactions of 
the sudden increase in temperature, to the extent shown by the 
curves of Figures | and 2. This holds up to the point where the 
two curves meet, after which the tubular reactor again becomes 
more advantageous. 

In the case of endothermal reactions, the CSTR is at a 
disadvantage from both the point of view of concentration and 
temperature. 
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Figure 3—Working curves for second-order endothermal 
reactions. 
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Figure 4—General comparison showing the variation of 
relative residence times with the conversion. 


A general comparison fot the 18 cases studied Is pre sented 
in Figure 4, obtained directly from the eighteen sets of working 
curves. It shows the relative residence times required to obtain 
given conversions. The shape of the curves, specially of those 
for exothermal reactions going through a minimum, is evident 


from consideration of the working curves as in Figures 1, 2 and 3. 


It is obvious from Figure 4 that within limitations the CSTR 
can be of definite advantage over the tubular reactor for exo 
thermal reactions. In isothermal reactions, and more so in 
endothermal, the reverse is true as the tubular can be decidedly 
preferable to the CSTR. 


Observation of the peculiar behavior for exothermal reac 
tions, where the CSTR is at times superior and then inferior 
to the tubular, has suggested that a combination of the two might 
give even superior results, taking best advantage of each re 
actor’s performanc¢ It has led to the study of the overall 


effects of combined reactors as outlined hereafter. 


COMBINED REACTORS 


Combinations of tubular and well-mixed reactors have been 
studied in the case of isothermal conditions’ for models 7M 
and M7 shown in Figure §. Noting that the effluent of the first 
section becomes the feed of the second, the overall effect for 
combined reactors can be determined tor adiabatic conditions in 
a way similar to that for isothermal, cither with respect to the 
conversion obtained for a given residence time or to the relative 


residence tme required for a given conversion 


(a) TM Reactors 


[he pertinent equations for cach section are presented in 
I abl 3. for first- and econd-orde! reaction I hose equation 
cannot be resolved into simple expressions as was the case for 
isothermal operation®® because of the compl x nature of the 
integral funcuons. Rather, a trial and error method is involved, 
based on the equations of Table 3, but making use also of the 


orking curves already mentioned 
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Figure 5—Illustration of model TM and MT reactors 
respectively. 


(1) Variation of conversion with mixing level at a given 

residence time. 

For a given residence time 6 or a total reactor volume such 
that the conversion be x, in a tubular reactor, the variation of 
conversion with the degree of mixing for any particular combina- 
tion TM is evaluated as follows. At a given value of m, the 
conversion x,’ at the outlet of the tubular section, as shown in 
Figure 5, is read directly on the proper working curve at a 
residence time 6,’ corresponding to that in the tubular section 
and equal to (1 m)@. Having thus taken care of the plug flow 
section one goes on to the particular expression of Table 3 for 
the CSTR section, noting that the term within brackets is 
obtained directly from the curve as it corresponds to 6, if all 
of conversion x took place ina CSTR. Using the proper CSTR 
working curve the procedure consists then in finding by trial 
and error the value of x such that the corresponding value of 
6, or C4,6 for a second-order reaction, will verify the expression 
for the CSTR section at the given value of # considered. An 
example 1s worked out in the Appendix. 


Of the many graphs prepared from such calculations, the 


one presented in Figure 6 1s an illustration typical of model TM. 

[he degree of conversion varies slightly at low mixing levels 
but may, depending on conditions, vary rapidly afterward as 
values of m increase. As expected, when R, increases, curves 
for second-order reactions approach more and more the curve 
for first-order. 

In the particular example given all conversions at complete 
mixing happen to be greater than that for m 0, in a tubular 
reactor. However, this is not the rule as lower conversions can 
be achieved in a CSTR; all depends, as can be seen in Figure 2, 
on the region where one happens to make comparisons. 


TABLE 3 


RELATIONS APPLICABLE TO MOvEL TM 


Order of 


Reaction 


Plug Flow Section CSTR Section 


| , (1 m)0 
Ca dx Ay, me 
o k(t v) ; \ ; Vy 
Rus 1 \ \ 
Ny 
Out 1 
x 
4 ( 190 1—m)¢ 1,0 C4,9u m ¢ 1,0 
xs dx \ vy,’ 
1 
yp RAI v) (CR; \ Ry. v) CR; ) \ 
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Figure 6—Variation of the conversion with the level of 
mixing for a TM reactor, under exothermal conditions. 


(2) Variation of relative residence times with mixing 
level at a given conversion. 


A more straightforward procedure than the one just described 
is available for this case. 

At a given x, the values of 6, and @, are read directly on the 
appropriate curves. A value of x,’ is then chosen determining 
on the tubular curve a corresponding value of 6,’.. From one of 
the expressions in Table 3 for the CSTR section, the value of 
6’ is then calculated from: 


The sum of 6,’ and @y’ gives the total residence time @ and m 
is thus easily obtained since m = 0y'/@. As to the relative 
residence time it is simply given by the ratio 6/6,. 

Two families of curves in Figure 7 show the variation of 
this ratio with the level of mixing. While very little difference 
exists in some cases between tubular, CSTR and combinations 
of reactors, there are others where it is rather pronounced. 
Here again, this factor depends on the region chosen for purposes 
of comparison. 


(b) MT Reactors 


For such reactors, as shown in Figure 5, the design equations 
work out to those presented in Table 4. 


Although methods similar to those used for 7Mf reactors 
could be employed in this case, a more convenient graphical 
method already described for isothermal reactions"® offers the 
advantage of being also applicable for exothermal and endother 
mal conditions. ‘The equations of ‘Table + are not even needed 
in this case but only the initial working curves, as shown here 

(1) Variation of conversion with mixing level at a given 

residence time. 


At a given value of m, the conversion x for a residence time 
# is obtained by following first the CSTR curve, as shown in 


PABLE 4 


RELATIONS APPLICABLE TO Mopet M1 


Order of 


R ; CSTR Section 
eaction 


Plug Flow Section 


’ vu | 4 fx dx 
Oy’ = m0 = (9, (l—m)6 
kay(l var’) | v, Ril v) 
i 
: Cay Ou’ = C4, : |¢ ‘Api C4,(1—m)0 
xu’ 1 dx 





ku( 1 —xy') (Re —xXu') vy! Rl v) (Rp— x) 
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Figure 7—Variation of the relative residence time with the 
level of mixing in a TM reactor, for two different conver- 
sions, under exothermal conditions. 
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Figure 8—Graphical construction to determine directly the 
conversion at a given level of mixing in an MT reactor. It is 
applicable to endothermal, isothermal and exothermal 
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Figure 9—Variation of the conversion with the level of 
mixing in an MT reactor, for endoethermal conditions. 


higure 8, up to point D. The latter, determined by @, . gives 
a conversion vy‘ and is easily obtained since @, me. A path 
DB is then followed, parallel to the portion Cad of the tubular 
curve, such that the abscissa of powut B corresponds to the 


Dp 


residence time 0 of the system. he ordinate ot point B deter 
mines the overall conversion tor the parucular combination 
considered 

In the above example the CSTR curve in Figure 8 lres below 


that tor the tubular. It may be noted however that the same 


195 





SECOND i 


/ 


/ 


MODEL MT 
X= 0.75 
EX OTHERMAL 


—— FIRST ORDER 


| | 





0.70 o75 080 085s 090 Oa 10 


Figure 10—Variation of the conversion with the level of 
mixing in an MT reactor, for exothermal conditions. 
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Figure 12—Variation of the relative residence time with 
the level of mixing at two different conversions, for exo- 
thermal reactions. 


construction is also applicable even if the CSTR curve happens 
to lie above the tubular one as in exothermal reactions. The 
method always consists in following the CSTR curve up to 6’ 
and then a path parallel to the tubular curve up to 6. 

Curves obtained by this method and showing the variation 
of conversion with mixing level are presented in Figures 9 and 
10 for endothermal and exothermal reactions respectively. 

In the endothermal case, the degree of conversion is seen to 
decrease gradually as the extent of mixing increases. The 
tubular reactor gives the highest conversion in all cases and the 
CSTR the lowest. 

For exothermal reactions, however, Figure 10 shows that a 
much higher conversion can be obtained at an intermediate mixing 
level than in either the tubular or CST reactors. 

The particular case of a second-order reaction, for R, = 1.25, 
is specially interesting as a conversion of 0.75 would be obtained 
in either a tubular or a CST reactor. Yet, a conversion approxi- 
mately 15% superior can be achieved at a mixing level of 0.6. 

The case of a first-order reaction is still more interesting as 
nearly complete conversion could be achieved at intermediate 
mixing levels while a tubular reactor would give one of only 
0.75 and a CSTR approximately 0.9. 

(2) Variation of Relative Residence times with mixing 

level at a given conversion. 

The graphical construction in Figure 11 is used to determine 
directly the residence time required to reach a given conversion 
at different mixing levels. 
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Figure 11—General graphical construction to determine the 
residence time required to reach a given conversion, at a 
certain level of mixing. 
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Figure 13—Comparison of the conversions obtained in MT 

and TM reactors at various levels of mixing. Full, dashed 

and dotted lines are for exothermal, endothermal and iso- 
thermal conditions respectively. 


Taking a given value of Ox’ determines a point /) on the 
CSTR curve at a conversion xy’. Following a path DA parallel 
to CB on the tubular curve gives at conversion x a point A 
which corresponds to the overall residence time 6. The latter 
together with Oy’ at point D yield the particular value of m. 
Since point B determines 6, for the conversion considered, one 

can obtain easily the values of the ratio 6/4, at all mixing levels. 

A plot of such values is shown in Fi igure 12 for two different 
conversions, in the case of exothermal reactions. 

No general conclusions can be drawn as to the best mixing 
level or MT combination to be used in such cases because the 
minima observed in the curve for the higher conversion are not 
found in those for the lower conversion. In this case a CSTR 
would be superior to any other arrangement. At higher conver- 
sion, however, a number of MT combinations would definitely 
provide a better utilization over a CSTR or a tubular reactor. 

The graphical method outlined for use with MT reactors is 
applicable because of the integrated function describing the 
operation of the tubular section and because the values of k 
and x corresponding to the lower limit of the integral are also 
those found in the effluent of the well-mixed section. 


COMPARISON OF TM AND MT REACTORS 

\ better picture of the relative behavior of the combinations 
7M and MT is perhaps obtained through a comparison of the two 
models for endothermal, isothermal and exothermal conditions. 

For illustration purposes, Figure 13 shows such a comparison, 
for a second-order reaction, giving the variation of m with x. It 
may be seen that for endothermal and isothermal conditions, 
while the conversion decreases for both models with increasing 
level of mixing, the 7M combination always gives a_ better 
performance at a given level of mixing than does the MT. 
Although no illustration is presented here as evidence, it may 
be mentioned that the horizontal spread between the 7M and 
WUT curves gets larger and larger as curves are drawn for 
conditions corresponding to increasing values of the conversion 
in a tubular reactor. 
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Figure 14—Comparison of models MT and TM showing, at 

a given conversion, the variation of relative residence times 

at various levels of mixing. Full, dashed and dotted lines 

are for exothermal, endothermal and isothermal conditions 
respectively. 


As to the exothermal case, model MT stands out as giving 
the optimum performance. On the contrary, far from improv ing 
performance, model TM gives poorer results. In the particular 
example presented, all MT combinations happen to be better 
than either the tubular or CSTR alone, and all the 7/1 poorer, 
because the conversions are the same at #7 = O and at m# = 1. 
General conclusions should not be drawn on that account but 
the behavior of model MT would be worthwhile investigating 
in any particular case at intermediate levels of mixing. 

Another basis of comparison, showing the variation of relative 
residence time with mixing level, is presented i in Figure 14. As 
expected, the same advantage of the MT combination is again 
evidenced for exothermal reactions. 

Since the curves of Figure 14 are for a large value of Ry, it 
may be noted that both curves for A/T and 7M combinations 
are identical for isothermal conditions as in the case of first-order 
reactions ®, 

A further method of comparing results would consist in 
plotting the conversion as a function of residence time for various 
levels of mixing. This has been done in Figure 15 for an iso- 
thermal second-order reaction. At a given level of mixing the 
curve for the 7M reactor lies closer to the tubular curve than 
does the one for the MT combination. It may be observed that 
as the conversion reaches a certain value, an increasing level of 
mixing would bring about a disproportionately larger ‘residence 
time for isothermal conditions, and more so for endothermal. 
The reverse would be observed on a graph presenting the 
exothermal case. 


Conclusions 


It has been shown that in exothermal reactions an optimum 
performance may be obtained with combined reactors of the 
MT model. Depending on the prevailing conditions a higher 
conversion would be obtained, or a smaller residence time 
required, at intermediate mixing levels. 

The usefulness of combined A/T reactors has already been 
mentioned in the course of investigations dealing with combus- 
tion, usual chemical reactions“®), the determination of mixing 
patterns® and incomplete mixing in a single reactor® as would 
be the case for instance in a tubular reactor at low flow rates‘2”). 

For the particular case of MT reactors, a graphical method 
is described which permits the direct and easy evaluation of the 
reactor performance. It constitutes in itself a valuable asset 
because it applies precisely to the combination of reactors which 
is susceptible of extensive applications. 

Although the present study has been limited in scope by the 
particular values of the kinetic parameters chosen, the interesting 
behaviour of MT reactors should lead to a more thorough 
investigation of the field which, backed by e xperimenti il evidence, 
would contribute further to the knowledge of reactor design. 
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Figure 15—Curves showing the variation of conversion with 
residence time. They permit the comparison of MT and 
TM reactors at different levels of mixing. 





Nomenclature 

A = refers to reactant A 

b = maximum temperature variation for complete conver- 
sion, °R. 

B- = refers to reactant B 

C = volumetric concentration, Ibs. /ft.* 

E = activation energy, B.t.u./Ib. mole 

k = reaction rate constant 

ky = reaction rate constant in the well-mixed section of com- 
bined reactors 

k,; = reaction rate constant in the tubular section of combined 
reactors 

m = level of mixing: fraction of the total volume of a combined 
TM or MT reactor which is perfectly mixed 

q = volumetric rate of flow, ft.*/hr. 

R- = universal gas constant, B.t.u./(Ib. mole) (°R.) 

Re = Ca,/Cap 

T = absolute temperature, °R. 

V = volume, ft.’ 

. Cy, = Ca 
x = conversion = ——— 
Cap 

xy = theoretical conversion ina CSTR 

x, = theoretical conversion in a tubular reactor 

vy’ = conversion in the CSTR section of a combined MT or 
TM reactor 

x;/ = conversion in the tubular section of a combined MT or 


TM reactor 
z = a constant 
6 = V/q = residence time, hr. 


Oy = residence time fora CSTR 
= residence time for a tubular reactor 
6,’ = residence time for the CSTR section in a combined MT 
or TM reactor 
0,/ = residence time for the tubular section in a combined MT 


or TM reactor 


Subscript 
i = refers to species 7 
* = refers to feed conditions 
M = refers to perfect mixing 
o = refers to initial conditions 
t = refers to plug flow operation 
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APPENDIX 


The following examples are worked out to illustrate the 
applications of the trial and error and graphical methods des- 
cribed in the text. They are based on the numerical values of 
the kinematic parameters presented in Table 1. 

(1) Find the conversion x, for a 7M reactor, given: exothermal 

first-order reaction, @ = 34.0 and m = 0.3. 

For m = 0.3, 6,) = 0.7 X 34.0 = 23.8. On the tubular 
curve of Figure 1, at a residence time of 23.8, the conver- 


sion x;) = 0.39 is read directly. By trial and error, a 

value of x = 0.79 is finally chosen on the CSTR curve 

at a res'dence time 6y = 20.1. This latter value when 
0.39 


iplied by (1 — x,’/x) gives Oy’ = 20.14 1 —- = 
multipl y ( /x) g Oy 079 
10.2 which checks the value of 64’ computed by the 
separate values of m@ = 0.3 (23.8 + 10.2) = 10.2. For 
the value of 6 = 34.0, a single tubular reactor would 
give a conversion of 0.75, while in this 7M reactor, one 
obtains a conversion of 0.79 at m = 0.3. 
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(2) Find the ratio 6/0, and the corresponding value of m for a 
conversion of 0.50 ina 7M reactor given: exothermal second- 
order reaction, Re = 1.25, and residence time in the tubular 
section corresponding to C48) = 10. 


On the tubular curve of Figure 2 at C40 = 10, x,’ = 0.145. 
At x = 0.5, on the CSTR curve, C4 Ou = 19.9. This 


x, 0.145 
value, when multiplied by { 1 — = ,orfi — “s) 


x 0.50 
gives C4,O0" = 19.9 & 0.71 = 14.1. The overall C40 = 
Oy’ 14.1 
14.1 + 10 = 24.1 and the value of m = —— on om O58, 
6 24.1 


As for x 0.50, C49 = 26.7 is read directly on the 


. 24.1 
tubular curve, the ratio 6/6, gives 967 = 0.90 for a value 


of m = 0.58. 


(3) Find the conversion x for an MT Reactor, given, endothermal 
second-order reaction, Re = 5, C40 = 12and m = 0.7. 

On the CSTR curve of Figure 3, a point D at xy’ = 0.67 
is located for C404’ = 0.7 X 12 = 8.4. Following the 
procedure outlined in Figure 8, the ordinate of point B 
gives the overall conversion x = 0.81. A single tubular 
reactor at C40 = 12 would give a corresponding con- 
version of 0.90. 

(4) Find the relative residence time and the level of mixing m in 
an MT reactor for a conversion of 0.90, given: exothermal 
second-order reaction, Re = 1.25 and C4 Om’ = 20. 

On the CSTR curve of Figure 2, a point D is located 
corresponding to C4 Om’ = 20 as indicated in Figure 11. 
A path DA parallel to segment CB of the tubular curve 
is then followed till the conversion of 0.90 is reached. 
The abscissa of point A reads C40 = 40.8. The value of m 
is easily computed as it equals 20/40.8 = 0.49. For a con- 
version of 0.90 a value of CaO = 47.5 is read directly on 
40.8 


47.5 





the tubular curve. The ratio 0/0, then becomes 


0.86 at this particular value of m = 0.49. 
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Speed of Sound in a Non-Ideal Gas’ 


WILLIAM SQUIRE? 


A method is proposed for computing the speed of 
sound in a non-ideal gas based on the approximation 


OH\ _H 
OE)” E 


s 


In conjunction with the usual generalized compressi- 
bility plot and an analogous reduced enthalpy correc- 
tion plot, it can be used near the critical point of a 
non-dissociating gas. The method is also applicable to 
highly dissociated and ionized gases. 


na recent paper by Buthod and Tien” a method was described 

for computing the speed of sound in a gas which deviated 
appreciably from the ideal gas law because of intermolecular 
forces. Their method is based on the use of generalized plots 
of the compressibility factor, Z, and its partial derivatives with 
respect to temperature and pressure. In this paper, an alternate 
approach is described which leads to a very simple approximate 
expression. This approximation is applicable both to conditions 
where the deviation from ideal gas behavior is due to inter- 
molecular forces and to conditions where the deviation from 
ideal gas behavior is due to dissociation and ionization. The 
second set of conditions is of importance in hypersonic flow 
investigations. 


Thermodynamic Basis 
The basis of the present approach is the thermodynamic 
identity 
oP P {0H 
Op }, op OE 


This identity is easily verified by using the well-known Bridgman 
table. It is obvious that 


OH H  _/OH/E 
eee BE ee Bi svete ova yf) 


dE), E dE 


and that the second term is a relatively small correction since 
the ratio of H/E is a slowly varying quantity. Therefore, the 
error involved in replacing (dH Ok), by H/E is not large and 
the convenient approximation 


‘Manuscript received October 12, 1960; accepted July 19, 1961. 

“Department of Applied Mechanics, Southwest Research Institute, San 
Antonio, Texas, U.S.A. Present address: Department of Aero-Space 
Engineering, West Virginia University, Morgantown, West Virginia, U.S.A. 
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PH : H 
uZ2 = a? a? ZRT —————_..... (3) 
pE H — ZRI 
can be used when adequate data along an isentrope is not available. 
However, we note that while the value of (0H/d0E), is 
independent of the reference state, the value of H/E is dependent 
on the choice of a reference state.* To obtain the correct answer, 
it is necessary to select the reference state so that H and E will 
vanish simultaneously. Fortunately, many tables of thermo- 
dynamic data are based on such a reference state. 


For a non-dissociating gas whose non-ideality is due to 
intermolecular forces, the usual corresponding state plots for 
Z as a function of the reduced pressure for various reduced 
temperatures can be used. A similar plot of a reduced enthalpy 
deviation, AH/T,., has been published by Edmister and can be 
used to evaluate H in Equation (3). This plot is reproduced as 
Figure 1 for the readers’ convenience. 


A gas such as air at a high temperature can be considered a 
mixture of ideal gases obeying the perfect gas law provided the 
correct molecular weight is used. The compressibility factor Z 

can be interpreted as the ratio of the nominal to actual molecular 
weight. This and the enthalpy are obtained from tables based 
on a solution of the equations for the equilibrium composition. 
Since the more recent tables“:>:" for air were prepared because 
of the growing interest in hypersonic flight, they include the 
speed of sound computed by a very accurate method. However, 
before they were av ailable, the author found E quation. (3) 
useful. It should also be useful for calculations involving hot 
combustion gases and other high temperature systems where 
complete tables are not available. 


Application to Ethylene 


Equation (3), in conjunction with the generalized plot for 
AH/T. (Figure 1) and a conventional Z plot, was used to 
calculate the speed of sound in ethylene at 7.9°K. (the critical 
temperature) as a function of pressure. In Figure 2, these 
values are compared with those calculated by Buthod and Tien 
and the experimental results cited in their paper. It is seen that 
the present method gives better agreement with experiment. 
The minimum in the vicinity of the critical point is reproduced 
though it is not as sharp as the experimental results indicate it 
should be. 


Apparently the error due to the approximation involved in 
Equation (3) is less than the improvement resulting from the 


*The author is indebted to Professor Agosta of Brooklyn Polytechnic 
Institute for an interesting correspondence on this point. 
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VALUES OF AH/T, , CAL/GM MOLE —°K 


G mame 
6 8S & S52 
REDUCED PRESSURE, P, 


Figure 1—Generalized enthalpy correction plot. 


TABLE 1 


ERROR OF THE APPROXIMATE FORMULA 
FOR SOUND SPEED OF HIGH TEMPERATURE AIR 


P = 10-% atm. 
(T in Exact value of % 
°K. X 10%) | Speed of Sound | error 
(meters/sec. ) 


P = 10° atm. | 

Exact value of | % 
Speed of Sound | error 
(meters/sec. ) 


2039 : 1517 
2207 : 1674 
2352 g 1816 
2684 ; 1973 
3209 ; 2163 
3700 i 2389 











use of the generalized enthalpy correction which has an accuracy 
comparable to Z instead of the derivatives of Z which are subject 
to a larger error than Z. Also, our formulation does not require 
the specific heat which is difficult to evaluate in this region. 


Application to Air at Very High Temperatures 


Table 1 compares the approximate speed of sound obtained 
by applying Equation (3) to the thermodynamic data for air 
between 6,000 and 12,000°K. compiled by Predvoditeley with 
the exact values. At 10~* atm., the air is almost completely 
dissociated at 6,000°K. and highly ionized, (43% electrons at 
12,000°K.). At 10% atm., there is an appreciable NO concen- 
tration and some oxygen dissociation at 6,000°K. and fairly 
complete dissociation, but only a little ionization at 12,000°K. 
It can be seen that over these widely varying conditions the 
approximate formula is frequently quite accurate, though there 
are regions where the error exceeds 5%. On the basis of a more 
extensive tabulation by Hochstim and calculations by Logan 
and Treanor“, it appears that the large error is associated with 
regions where the composition is changing rapidly rather than 
with a high degree of dissociation or ionization. 


PRESENT THEORY 


(approaches @ between 
1.14 @ 1.16) 


° 
° 


° 


SPEED OF SOUND,METERS/ SEC 


© EXPERIMENT (HERGET) 


0.20 0.40 0.60 0.80 


PRESSURE 
CRITICAL PRESSURE 


Figure 2—Variation of the speed of sound in ethylene, at 
the critical temperature, as a function of pressure. 
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Nomenclature 


a = conversion factor to express u, directly in meters/second 


, meters gram \'? 
It is equal to 31.83 - 
second \liter-atm 
= specific heat at constant pressure (calories/grams °K.) 
specific heat at constant volume (calories/grams °K.) 
internal energy per unit mass (calories/grams ) 
enthalpy per unit mass (calories/grams) 
pressure (atmospheres) 
gas constant (liter atmospheres/grams °K.) 
enthropy per unit mass (calories/grams °K.) 
temperature (°K.) 
speed of sound (meters/second) 
= compressibility factor (dimensionless) 


= reduced enthalpy correction, enthalpy difference divided 
rT. by critical temperature (calories/grams °K.) 
p density (grams/liter) 
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The Cai 


Kinetics of the Vapor-Phase Oxidation of 
Toluene over a Vanadium Catalyst’ 


J. DOWNIE’, K. A. 


SHELSTAD® and 


W. F. GRAYDON‘ 


Rates of oxidation of toluene were measured as 
a function of temperature and reactant concentration. 
A fixed bed reactor was operated at low conversions 
over the following range of conditions, 300°C. to 
350°C., 1 x 10° to 15 x 10° moles oxygen per litre 
of gas, and 1 x 10° to 6 x 10° moles toluene per litre 
of gas. 


The rate data for the oxidation of toluene were 
correlated by the following equation, which is based 
on a model by Hinshelwood. 


ky k, C; c. 
n= 
. ee, & FAS. 


The values of k, obtained agree with those re- 
ported previously for naphthalene oxidation"). 

It is concluded that the Hinshelwood model is a 
useful one for this type of reaction, and within the 
limits of available data provides a satisfactory inter- 
pretation. 


yo the past decade interest has developed in vapor-phase 
catalytic oxidation which has been of a more fundamental 
nature than had hitherto been evident. E pirical and engineering 
advances had placed many processes on a sound operating basis, 
and, in doing so, had outdistanced the understanding which 
results from more fundamental approaches. 


Several basic studies have now been made, especially with 
regard to the catalytic oxidation of naphthalene. These have 
been discussed in a previous paper“. However, data have not 
yet been presented for the dependence of the rate of oxidation 
of toluene on reactant concentrations and temperature. 


Experimental Work 


The apparatus, shown as a flow diagram in Figure 1, con- 
sisted of independent flow systems, for oxygen and nitrogen, 
leading to a mixing chamber from which were drawn primary 
and secondary supplies, Toluene saturation of the primary 
supply and reunion with the secondary supply preceded the 
reactor. The reactor was kept at constant temperature and the 
gases passed through a thin catalyst bed. The products were 
collected in suitably designed liquid air traps. ‘The apparatus 
was operated at low conversions to study the rate of oxidation 
directly as a function of known inlet concentrations. 


One sample of the reaction products was analysed spectro- 
photometrically, another by potentiometric titration. From the 
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results the moles of p- benzoquinone, benzaldehyde, benzoic acid 
and maleic anhydride in the reaction products were calculated. 
This was the information which was used to calculate the moles 
of toluene oxidized per run, and hence the rate of oxidation of 
toluene. 

Carbon dioxide was measured by the weight increase of a car- 
boxite trap. After finding negligible amounts of complete oxida- 
tion the carbon dioxide checks were only used intermittently. 

From the product distribution of the low conversion runs it 
may be seen that of we ne oxidized approximately 86% 
reacts to benzaldehyde, 7% reacts to p-benzoquinone, 3% reacts 
to maleic anhy dride, and a reacts to benzoic acid. The error 
in the estimated rate of oxidation, due to analytical errors, is 
therefore approximately that of the benzaldehyde determination, 
i.e., 4% 


Theory 

The reaction takes place at the solid surface, but our informa- 
tion about reactant concentrations refers to the gas phase. 
Therefore we require a model which permits the calculation 
of surface concentrations in terms of gas phase concentrations. 

In the past many reactions have been explained by assuming 
adsorption equilibrium to be established. Then L angmuir’s 
theory relating known concentrations to surface concentrations 
has been used. This approach has been most successful for high 
temperature reactions, ¢.g., incandescent metal wire catalysts, 
but there is another class of reactions where the reaction temper- 
ature is relatively low. For the latter Taylor found a slow rate 
of adsorption®. When the rate of adsorption of reactants is 
of the same order of magnitude as the rate of chemical reaction 
the assumption of adsorption equilibrium is no longer valid, and 
rate relationships based on Langmuir-Hinshelwood or Rideal 
mechanisms are on a doubtful footing. 

The present study involved using a vanadium oxide catalyst 
at 300° to 350°C., and it seems reasonable to assume adsorption 
equilibrium was not established. 

It was decided to use the simplest steady state situation 
first and to check the resulting relationship against the experi- 
mental data. 


COLD TRAP 


REACTOR 


MIXING 


CHAMBER 
CARBURETTORS IN 


A CONSTANT 
TEMPERATURE 
BATH 


FLOWMETER 


Figure 1—Flow diagram of the apparatus. 





This treatment follows Hinshelwood’s general treatment on 
page 208 of “Kinetics of Chemical Change”. The introduction 
of N is for the general reaction of this particular study. 


The following assumptions are made, 


1. Only oxygen is adsorbed on the catalyst surface. 
2. Reaction occurs when toluene in the gas strikes adsorbed 
oxygen 
. The rate of desorption of oxygen is negligible. 
. A steady state 1s established in which, Rate of removal 
of oxygen by chemical reaction = Rate of adsorption of 
oxygen. 
C73 + (N) O2 — Products 
Therefore, 
Rate of removal of oxygen by Chemical reaction 
= N X Rate of oxidation of toluene 
NR&,C,S 
Therefore, 


NRCS =k. C.(1 — S) 


bt, 
Lt Nn k, c 
Therefore, 


Rate of oxidation of toluene = 


Discussion of Results 

The rate of catalytic oxidation of toluene was determined as 
a function of toluene concentration, oxygen coricentration, and 
catalyst bed temperature. This involved measuring the rate of 
oxidation both in the presence and in the absence of catalyst. 
Knowing the oxidation in the absence of catalyst under ‘the 
conditions of the catalytic run, the oxidation due to the presence 
of the catalyst was found by difference. Details are given in 
reference 4. 

The experimental results presented in Figures 2, 3, 4, 5 and 
6 have been corrected for oxidation in the absence of catalyst, 
and have been brought to a constant concentration basis for the 
parameters. The solid lines were calculated from Equation (1). 
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INDUSTRIAL CATALYST 
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ow 
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D - (CyHg)* 1.0 X10 MOLES/LITRE 
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Du 


4 6 12 
MOLES OXYGEN PER LITRE OF GAS Xx 10> 


Figure 3—Rate of oxidation of toluene as a function of the 
concentration of oxygen, at 350°C., using the industrial 
catalyst. 
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LAB. CATALYST 
350°C 
(C7Hg)= 5 X 10°4 MOLES/LITRE 


Nn w 


MOLES TOLUENE OXIDIZED PER g. CATALYST PER SECOND xX 108 





% 4 8 12 


MOLES OXYGEN PER LITRE OF GAS X 103 
Figure 2—Rate of oxidation of toluene as a function of the 
concentration of oxygen, at 350°C., using the laboratory 
prepared catalyst. 


Equation (1) becomes suitable for use in correlation when 
inverted. 
1 1 


mn kG 


Equation (2) was used to test the effectiveness of the 
Hinshelwood interpretation in the correlation of the experimental 
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350°C 
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Figure 4—Rate of oxidation of toluene as a function of the 
concentration of toluene at 350°C., using the industrial 
catalyst. 
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Figure 5—Rate of oxidation of toluene as a function of the 
concentration of toluene, at 325°C., using the industrial 
catalyst. 


results. The straight line relationship was a good fit of the data. 
k, and k,/N were calculated from a least squares treatment of 
the data and are given in Table 1. 

The average values for k, and k,/N were used in Equation (1) 

to calculate the lines shown in Figures 3, 4, 5 and 6. The average 
relative errors of the calculated rates were estimated for each set 
of data and were as follows: 
Set No. 1 2 3 4 5 
Average % Error t 3 7 6 10 
Inspection of the curves shows a satisfactory correlation of the 
data, indicating the relative insensitivity to the small variations 
in k, shown in Table 1. N (the moles of oxygen used per mole 
of toluene oxidized) was determined, from product distribution 
data, to be 1.35 throughout the range of conditions in Sets | to 5. 
k, was then calculated and an Arrhenius-type plot used to obtain 
Arrhenius energies of activation for k, and k, (Figure 7). 

These results can be summarized as follows: 

the catalytic oxidation of toluene over industrial catalyst | 
in the temperature range 300° to 350°C. can be described by 
the following relationship: 


where N = 1.35 
k, = 1.39 X 10-4 at 350°C. 
5.22 X 107* at 325°C. 
1.78 X 10-5 at 300°C. 
Arrhenius activation 
29,400 cal./mole 


energy of k, = 
k, = 2.38 X 10-* at 350°C. 

1.27 * 10-4 at 325°C. 

3.76 X 10-5 at 300°C. 


Arrhenius activation 
26,400 cal./mole 


This treatment has been applied to the rate of oxidation data 
for naphthalene, and we are now in a position to show that the 
theory not only is useful in correlating the studies independently 
but provides a link between them. 


energy of k, 
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Figure 6—Rate of oxidation of toluene as a function of the 
concentration of oxygen, at 300°C., using the industrial 
catalyst. 


According to the above interpretation k&, should be independ- 
ent of the organic compound oxidized when the same catalyst is 
used. Since the same industrial catalyst was used in both toluene 
and naphthalene work, a check is available. 

At 350°C. 
k, from naphthalene = 0.93 X 1074 
k, from toluene = 1.39 < 10 


At 300°C. 


k, from naphthalene 
k, from toluene 


1.30 X 10-5 
1.78 X 10-5 


Arrhenius energy of activation of k, 


from naphthalene = 28,100 cal. mole 
from toluene = 29,400 cal. mole 


The &, values from toluene and naphthalene differed by a factor 
of 70. Therefore, the small difference in k, shown above was 
taken as supporting the view, arrived at on the basis of the main 
body of results presented previously, that the Hinshelwood 
approach is a useful one for this type of reaction and within the 
limits of the available data provides a satisfactory interpretation. 

Mars and van Krevelen® arrived at similar conclusions for 
their work on benzene, naphthalene, and anthracene, although 
they state that they could not include toluene. They did not 
publish results on toluene and therefore no comparison can be 
made. 

TABLE 1 
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Figure 7—Arrhenius-type plot of the specific rate constants 
as functions of temperature. 


The results of a series carried out using a laboratory prepared 
catalyst are shown in Figure 2. The results were of the same 
form as those using the industrial catalyst. 


A series of oxidation runs was started to study the effect of 


contact time on product distribution. The results are shown in 
Figure 8. 

The results show a levelling off in the yields of benzaldehyde 
and p-benzoquinone, and an increase in yields of benzoic acid 
and maleic anhydride over the entire range studied. The exit 
concentration of each of the above products was calculated and 
examined as a function of contact time. 

These results were used in the following speculative calcu- 
lations. 

The rate of oxidation of benzaldehyde was calculated assum- 
ing that all p-benzoquinone, maleic anhydride, and benzoic acid 
resulted from the catalytic oxidation of benzaldehyde. These 
calculations were carried out for the four runs at the two lowest 
contact times in order to minimize the error due to CO. forma- 
tion. The total conversion of toluene under these conditions 
was 4% and 9%. 

The rates thus calculated were used together with benzalde- 
hyde concentrations to determine k,. ‘The value found in this 
way was 

k, = 0.66 X 10-4 at 350°C. 
i.c., about half the value previously obtained from toluene 
oxidation data. 

This calculation was repeated for benzoic acid oxidation 
and resulted in a k, value of 0.27 & 10~4 which is abcut one-fifth 
of the toluene value. 

If, as is possible, there was some total oxidation during the 
9% conversion runs, accounting for this will result in k, values 
which approach the k, value from toluene oxidation data. 

Therefore, as a tentative conclusion, it may be possible to 
extend the treatment of the main body of results to cover 
benzaldehyde and benzoic acid oxidation 


Conclusions 

The rate of oxidation of toluene over a vanadium oxide, 
potassium sulphate promoted, catalyst, can be described by the 
following relationship: 
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Figure 8—% yield of reaction products as a function of 
W/F. 


RAC, C, 
OT ht + OA 

Values of k, from data on the oxidation of naphthalene” 
agree to within 33% of those from the present work, and values 
of the temperature caulbeiein of k, agree to within 4%. This 
interpretation therfore provides a link between oheme and 
naphthalene oxidation on the same catalyst. 

The Hinshelwood model is a useful one for this type of 
reaction and within the limits of available data provides a 
satisfactory interpretation. 

Further oxidation work on toluene at varying contact time 
indicated that it may be possible to extend this approach to 
benzaldehyde and benzoic acid oxidation. 
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Nomenclature 

C = gas concentration, g. moles per liter 

k, = specific rate of adsorption of oxygen 

k, = specific rate of reaction of toluene 

N = moles of oxygen used per mole of toluene oxidized 

r = reaction rate, moles of toluene oxidized per g. of catalyst 


per second 


S = fraction of surface covered by oxygen at steady state 
T = temperature 

Subscripts 

0 = oxygen 

t = toluene 
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Determination of Vapor-Liquid 


Equilibrium Data from Total Pressure- 
Liquid Composition Measurements’ 


JAMES C. K. HO*?, O. BOSHKO? and 
BENJAMIN C.-Y. LU? 


Digital computer methods are developed for 
evaluating equilibrium data from total pressure-liquid 
composition measurements of volatile binary non- 
electrolyte mixtures. The Redlich-Kister equation is 
employed for representing liquid activity coefficients. 
The constants B, C and D are evaluated by the method 
of least squares for non-linear formulas. Equilibrium 
vapor compositions are also evaluated by means of a 
step-wise integration procedure. A comparison of both 
methods is presented. 


+ ie pressure-liquid composition measurements are in general 
easier than the equilibrium vapor composition measure- 
ments, especially for solutions of high relative volatilities and 
for solutions which are difficult to analy ze. Several methods 
have been proposed in the literature for evaluating liquid 
activity coefficients, partial pressures and equilibrium vapor 
compositions from the total pressure-liquid composition measure- 
ments“:?-3-4..6.7), In general, these methods involve evaluations 
of either constants of equations representing the liquid activity 
coefficients or equilibrium vapor compositions directly by means 
of a step-wise integration procedure. All methods involve, 
however, the Gibbs- ‘Duhem equation. 

In the course of inv estigation of various thermodynamic 
properties of non- -clectrolyte solutions in this laboratory, it 
appeared desirable to develop a method for evaluating excess 
thermodynamic functions which are represented by a power 
series with respect to liquid compositions. The use of a digital 
computer (IBM 650) was intended. In this paper, application 
of this method is illustrated for evaluating activity coefficients 
from total pressure-liquid composition measurements. It is also 
of interest to this laboratory to evaluate equilibrium vapor 
compositions directly from the total pressure measurements 
without involving any conversion to liquid activity coe ficients. 
The advantages of these two approaches are compared and 
discussed in this paper. 


Evaluation of Activity Coefficients 


Total pressure-liquid composition measurements are usually 
made at constant temperature conditions. The Gibbs-Duhem 
equation of the form 


0 log y1 0 log 2 
x — + Xe - =(@... (1) 
’ Ox, oy Ox, . 


is generally employed for developing equations for representing 
liquid activity coefficients. It is obvious that Equation (1) is 
applicable only at constant temperature and pressure and all the 
equations developed from it are applicable only at these con- 
Ghee cecuscncepneceendenncsneadennsasendenneucnncscenancecee cescncsecasoncncaes 
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ditions. However, it is due to the fact that the pressure effect 
on liquid activity coefficient is expressed by 

0 log ¥: 0; — 0; Q2 
ap“ papaey eee OD 


which is small and can be neglected over a moderate range of 
pressure, that integrated forms of the Gibbs-Duhem equation 
can be employed at isothermal conditions. 

The expression for representing the excess free energy 
function, Q, as proposed by Redlich and Kister‘® is employed 
in this investigation 


It may be readily seen that from Equations (3) and (1) the 
i activity coefficients may be represented by the following 
expressions: 


logy: = (1 —»1)*>[B + C(4x1 —1) + D(2m, —1)(6x, -1)+...)]. (4a) 
and 
log y2 = x7 {[B + C(4x,—3) + D(2x1—1) (6x,—5) + ...]..(4b) 
If the vapor phase is considered ideal, then 
Pe i BE So iavcsstante (5) 


Simultaneous solution of Equations (4) and (5) would provide 
the required information. The problem is therefore to determine 
the constants involved in Equations (4a) and (4b) so that the 
total pressure values calculated from Equation (5) agree with 
the experimentally observed values. 

Christian employed expressions similar to Equations (3) 
and (4) for representing liquid activity coefficients and deter- 
mined the constants (restricted to three) using three experimental 
points. The three points were chosen either from the smoothed 
curve or from three actual data points without any preference 
in the choice of liquid compositions. Should only three constants 
be considered in the equations representing activity coefficients 
and only three points be employed, it was observed in this 
laboratory‘ (% that the - three points would be at x; = 0.25, 
x, = 0.5, and x, = 0.75. The evaluation of the constants could 
be very much simplified. Such a computing program was 
developed in this laboratory”. 

In the present investigation, the method developed uses all 
the experimental points in the evaluation of the constants. The 
principles involved could be employed for any number of con- 
stants involved in Equations (4a) and (4b). The IBM program 
developed for this paper is limited however to three constants, 
B, C and D. The mathematical approach is similar to that of 
Baker™ and outlined below. 

As the Equations (4a) and (4b) are not linear, the constants 
B, C and D cannot be determined directly by means of the method 
of least squares. A set of initial values of B, C and D may be 
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assigned, using constants of a similar system as a guide or 
approximated by assuming that the solution behaves as a regular 
solution so that C and D are zero. Using the first approxi- 
mation, the pressure residuals, AP, and the derivatives 0P/OB, 
@P/8C and OP/AD are determined by means of the following 
equations. 





oP OP = 
= evd. — Praca. = AB — + AC— + AD—...... 6 
AP = Povsva. led apt aC (6) 
ap 
— = 2.303 [xi(1 = x)? pr "1 - x? (1 — xX) po ae ai we a3 
OB 
oP 
5c 7 2303 [a(1 — a)* 4s — 2 pn + 
x? (1 — x1) (4x; — 3) p° 72! Cee eoeseseeseecsees (8) 
and 
oP 
aD = 2.303 [x, (1 — 1)? (2x. — 1) (6m. — 1) poy + 
x1? (1 — 2%) (2x; — 1) (6x; -_ 5) p° Val ss os aes .(9) 


The activity coefficients y, and 2 in these expressions are 
represented by Equations (4a) and (4b). The procedure of 
calculation is to determine the changes AB, AC and AD in B, 
C and D, so that the pressure residuals be reduced to a minimum 
by means of the least square method, which is described below. 

These increments are added to the initial values of B, C and D 
for improved approximations. Thus, the quantities B, C and D 
are determined by successive approximations. The process is 
repeated until the increments AB, AC and AD are equal to or 
less than the tolerance desired in the values of B, C and D. 


For n number of observations, we have n residual equations: 


OP, aP, aP, 

Bo + Acs + ADS = ar, 
OP, OP, OP. 

AB? 4 Ac Ap 2 ~ ap, 
mn” “te * "SS : 
OP,, OP, P, 

- ; )- = AP, 

ABS" + ACS. + AD So A 


Therefore, the “normal” Equations‘ are 
P\? OP\ (AP 
an ¥(35) + 4¢2(35) (sc) + 
OP oP oP 
D 7 “g 
” =(55) ol Lal (55) tai 
OP\ {AP 

ae £(55) (3c) + 8¢ (3c) 

1s =(52) ( 
and 


oP OP 
B 
. =(55) (5 5) + ac ”s 
oP? oP 
D = Dy ee ie ei Se d 
. “(t) 2 a(t) - 


where the summations }~ are taken over all m observations. 
There are three normal equations as there are three quantities 
B, Cand D to be determined. Multiplying the Ist, 2nd, .... and 
nth residual equation by AP,, AP:, .... and AP, respectively, 
and adding all the corresponding factors together, Equation (13) 
is obtained, 


a" 


OP 
> AP = ae ; eee) 
(57) ( 
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AB (AP) ) + AC X(AP) (5) + 
AD >-(AP) (55) - RMS Biases eka aeens (13) 


from which one may judge the quality of the fittings. 

In the calculations, it is observed that the values of first 
approximation do not affect much the computing time required. 
For the system ethanol-water at 25°C.“ two sets of initial 
values of B, C and D were assumed. In the first set, B = 0.200, 
C = 0.040 and D = 0.020 were employed and in the second set 
B = 0.830, C = 0.000 and D = 0,000 were employed. The 
difference in the required machine time is less than half a minute. 

There are two modes of output in the program. In the 
condensed output, the following items are printed: 


(1) Initial estimates of the constants, B,, C, and D, and the 

tolerance desired. 

(2) Final values of B, C and D and number of iterations 

performed. 

(3) The quantities, x, Popsva., AP, pi, Po, Yi and Y2. 

In the optimum output, the quantities AB, AC, AD and AP are 
printed for each iteration. The course of the computation may 
be easily followed step by step. The program was originally 
written in FORTRAN symbolic language. The machine time 
required per set of data depends, in general, on the number of 
experimental data, initial estimates of B, C and D, and the desired 
tolerance. As an average, four iterations are adequate for a 
tolerance of the order of 10~* and approximately four minutes 
are required per set with condensed output. The program 
occupies 1,005 locations. 

The calculated B, C and D values for 13 systems at 20 
temperature conditions are listed in Table 1. A typical output 
is shown in Table 2 for the system ethanol-water at 25°C.“), 
The computed equilibrium vapor compositions of this system 
are compared with the experimental values in Table 3. 

It is of interest to note that in Table 1, with the exceptions of 
the data for the extremely imperfect systems ethanol-isooctane* 
and ethanol-methylcyclohexane, the |AP\,,. values are gener- 
ally small, indicating that the three constant Redlich-Kister 
equations are adequate for representing liquid activity coefh- 
cients. This is significant in testing data consistency. It 
may be mentioned that for the system ethanol-water at 74. 79°C, 
the computed B, C and D v alues from the total pressure measure- 
ments (Table 1) are in good agreement with the values calcu- 
lated“ from the experimentally determined vapor-liquid 
equilibrium compositions“®., However, in the case when 
experimental data were used, the values B = 0.5311, C = 
—0.1650 and D = 0.0441 were obtained after the data were 
adjusted for consistency. The advantage of evaluating activity 
coefficients using total pressure measurements is obvious. 


Evaluation of Equilibrium Vapor Compositions 


In the above section, the evaluations of the constants B, C 
and D of the Redlich-Kister equation, the activity coefficients 
and the partial pressures have been described. ‘The equilibrium 
vapor compositions may be determined by using Equation (5). 


Thus 

AZUANfPr /P = p/P........ 
and 

Yo = X22 p/P = pro/P = 1 —m.......... (14b) 


However, it is also of interest to develop a computing method 
for evaluating equilibrium vapor compositions directly from the 
total pressure measurements without involving any conversion 
to liquid activity coefficients. 

The relationship between the total pressure measurements 
and the compositions of liquid and vapor phases may be expressed 
by 

dP uM — MH 


= — dy,..... te 
F ym(1 - yn) (4 
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TABLE 1 
CALCULATED CONSTANTS AND DEVIATIONS 



































|AP| av. 
System Temp. °C. B Cc D mm. Hg. 

1. Acetone-Chloroform(® 25 —0.4276 0.1105 0.0571 0.71 
35 —0.3343 0.0701 0.0357 0.47 

50 —0.3877 0.0941 0.0436 0.33 

2. Benzene-n-Heptane® 80 0.1578 0.0455 0.0142 0.28 
3. Carbon Disulfide-Benzene®” 20 0.1623 0.0594 —0.0406 0.24 
4. Carbon Tetrachloride-Benzene® 40 0.0538 0.0000 —0.0024 0.03 
70 0.0469 —0.0002 0.0010 0.05 

5. Chloroform-Ethanol®» 44.98 0.4540 0.2295 0.0214 1.47 
6. Ethanol-Isooctane@* 25 0.9838 —0.0740 0.2270 1.12 
50 0.9458 —0.0839 0.2610 5.75 

7, Ethanol-Methylcyclohexane 35 0.9498 —0.1074 0.2358 1.87 
55 0.9185 —0.0867 0.2260 4.78 

8. Ethanol-Toluene 35 0.7946 —0.1149 0.1962 1.27 
55 0.7596 —0.0897 0.1188 1.53 

9. Ethanol-Water(” 25 0.5251 —0.1270 —0.0258 0.25 
rs , as) 74.79 0.5246 —0. 1595 0.0476 1.80 

10. n-Hexane-Chloroform 65 0.2500 —0.0150 —0.0101 0.67 
11. n-Hexane-Ethanol@® —10 1.0007 0.1194 0.1947 0.18 
12. Methanol-Water®” 25 0.2134 —0.0218 —0.0603 0.37 
13. Methylcyclohexane-Ethanol?® 0 1.0149 —0.1376 0.1291 0.33 





in which the vapor phase is assumed to be ideal. Equation (15) 
may be obtained by substituting the following expressions 


Aa =n P, dp, = Pdy, + dP 
and 


pe = i, dpe = Pdy, + 2 dP 


in the Gibbs-Duhem equation of the following form 


0 log pi 0 log br 
x ——— Xe = eae eae ume 16 
” ( Ox, )., : ( dx, T.P , ie 


Equilibrium vapor compositions may be evaluated from Equation 
(15) by means of a step-wise integration procedure as illustrated 
by Othmer et al®. In the integration, the smaller the increment 
of Ax used, the better are the integrated results. For a small 
increment of Ax, for example, 0.01, it would be a quite time- 
consuming task to do the calculation using a desk calculator. 
An attempt has been made in this laboratory to develop a satis- 
factory digital computing program for doing such a calculation. 





Without a computer, the step-wise integration procedure 
may be carried out by using Equation (15) in the rearranged 
form of 


sn ha. Mo US ees weQh) 
yi(1 —_ 41) 
TABLE 2 


COMPUTED VALUES FOR SYSTEM ETHANOL-WATER AT 25 °C. (12) 











No. of | 

B Cc D Iterations 
0.5251} —0.1270} —0.0258 + 

x P mm.Hg.}Pobsva — Peated pr pe v1 v2 
0.0523} 33.17 —0.34 10.89 22.62 | 3.528 | 1.005 
0.0917 38.44 —0.24 16.79 21.89 | 3.103 | 1.015 
0.1402 43.42 0.16 22.13 21.13 | 2.675 | 1.035 
0.1670 45.69 0.54 24.39 20.76 | 2.475 | 1.050 
0.2027 47.21 0.04 26.85 20.32 | 2.245 | 1.073 
0.2849 50.33 0.04 30.84 19.45 | 1.834 | 1.145 
0.3368 Sti] —0.42 32.63 18.96 | 1.642 | 1.204 
0.4902 53.95 —0.30 36.84 17.41 | 1.274 | 1.438 
0.5820} 55.71 0.18 39.47 16.06 | 1.149 | 1.618 
0.7810 58.08 0.24 47.13 10.71 | 1.023 | 2.058 
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and a plot of Jn P vs. x;. The integration starts at x; = 0, where 
y: = Oand P = p.°. In order to carry out the same procedure 
by means of a computer, it is necessary to represent the Jn P vs. x; 
plot analytically. However, as the representation of the function 
In P by means of a polynominal expressed in terms of x; was found 
to be not sensitive enough to the variations of x,, it was decided 
to represent the total pressure P by means of a polynominal 
expressed in terms of x;. Thus, 


P =a + bx, + cx? + dx? + ext........... (18) 


in which a = po° anda +b +¢ +d +e = p,°. In other 
words, for total pressure measurements obtained at isothermal 
conditions, Equation (18) involves three unknown constants. 
Should more constants be employed in Equation (18) the 
computing time required for performing the step-wise integration 
procedure would be more than the time required for ev aluating 
the activity coefficients and partial pressures. It was observed 
also that for azeotropic systems having rather flat total pressure 
readings in the vicinity ‘of the azeotropes, Equation (18) was 
not adequate for re presenting the data. The equilibrium composi- 
tion curves obtained by starting the integration from x, = 0 and 
x; = 1 respectively would not meet at the azeotrope. Although 
the integration procedure could be carried out initially at the 
azeotrope, the time required for locating the azeotropic composi- 
tion by solving dP/dx = O was considerable and the result 
obtained was misleading due to the poor fitting. It is because 
of these difficulties that the computing program developed is 
limited to non-azeotropic systems. For azeotropic systems, it 
is recommended that the method presented in the first section of 
this article be followed. 


TABLE 3 


COMPUTED AND EXPERIMENTAL EQUILIBRIUM VAPOR 
COMPOSITIONS } FOR THE SYSTEM ETHANOL -WaTER AT 25°C, 





v1 Yexpt. | Yealed. 








| 
0.0523 0.3166 0.3283 
0.0917 0.4334 0.4367 
0.1402 0.5129 0.5097 
0.1670 0.5450 0.5338 
0.2027 0.5687 0.5687 
0.2849 0.6106 0.6128 
0.3368 0.6287 0.6377 
0.4902 0.6791 0.6829 
0.5820 0.7096 0.7085 
0.7810 0.8161 | 0.8115 








0 0.2 0.4 0.6 0.8 1.0 
X| 


Figure 1—The calculated and experimental equilibrium 
composition values for the system methanol-water at 25°C. 
O Experimental; — Calculated. 


In the computing program, Equation (15) is rearranged in 

the following manner: 

AP be 

A i, ce capsird seis 

I n(1 — 1) 
An increment of x; = 0.01 is employed in the step-wise calcula- 
tions. It should be mentioned, however, that the computing 
program developed is not limited to this increment. Between 
x, = 0, y: = O and at the end of the first increment, 

AP = Ps, = 0,01 — Py, = 0 

Ay =H 

and 
AP/P +x 
1 + AP/P 





nee 


The equilibrium vapor composition at the end of the first 
increment is thus obtained. At the end of the second increment, 
the equilibrium vapor composition may be obtained by solving 
the quadratic equation for 9; resulting from E quation (19). The 
larger value of the two roots is taken to be the answer. This is 
necessary as jy; is always greater than x for non-azeotropic 
systems. Subsequent integrations are carried out in the same 
manner. As an average, the machine time required per set of 
data is approximately four minutes, which is about the same as 
the time required for evaluating the values of B, C and D. The 
integrated values of the equilibrium vapor compositions from 
the total pressure measurements for the system methanol-water 
at 25°C.“ are compared with the experimental values in 
Figure 1 as an illustration. 


Conclusions 


Two digital computer methods are developed for evaluating 
vapor- liquid equilibrium data from total pressure-liquid composi- 
tion measurements of volatile binary non-electrolyte mixtures. 
One of which calculates the constants B, C, D of the Redlich- 
Kister equation, hence the values of activity coefficients, partial 
pressures and equilibrium vapor compositions (Method 1). The 
other evaluates the equilibrium vapor compositions directly by 
means of an integration procedure (Method 2). Although 
Method 2 does not involve any conversion to liquid activity 
coefficients and offers the equilibrium vapor compositions 
directly, the machine time required per set of data is about 
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the same as that required by Method 1, should an increment of 
0.01 of x is chosen. Also, Method 2 developed is limited to 
non-azeotropic system. Method 1 offers, on the other hand, 
complete equilibrium values. 


The calculated results using Method 1 provide further 
support of the suggestion“® that three-constant Redlich-Kister 
equation is adequate for representing liquid activity coefficients 
for most systems. 


Acknowledgment 


The authors (J. C. K. Ho and B, C.-Y. Lu) are indebted to the 
National Research Council of Canada for financial support. The authors 
are also grateful to the Computing Centre of the University of Ottawa 
for making the computer time available and for its assistance in develop- 
ing the computing program. 


Nomenclature 
B,C,D = constants 
a,b,c,d,e = constants 
P = total pressure 
p = partial pressure 
p° = vapor pressure of pure component 
Q = excess free energy function 
R = gas constant 
7 = absolute temperature 
5 = partial molal volume 
v = molal volume 
x = liquid composition, mole fraction 
y = vapor composition, mole fraction 
Y = liquid activity coefficient 
Subscript 
1,2,7 = components 
1,2,...n = no. of observations 
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The Hydraulic Characteristics of 


Two-Phase Co-Current Flow through 
Glass Fibres’ 


D. J. GUNN? and A. R. AITKEN® 


The mechanism of flow of air and water through 
packed glass fibres has been studied. Pressure 
drops and liquid hold ups were measured as a func- 
tion of gas and liquid flow rates and the packing 
characteristics. It has been found that the pressure 
difference measured is dependent upon the history 
of previous gas and liquid flow rates in the bed. The 
liquid hold up is defined entirely by the gas and 
liquid flow rates. An analysis of the flow mechanism 
shows that the capillary forces interact with the fluid 
pressures in the bed. For example an increase in gas 
flow rate followed by a reduction to the same flow 
rate has been shown to reduce the resistance of the 
bed to gas flow even though the liquid flow rate is 
unchanged. A measure of the irreversibility of the 
flow operation is suggested. 


co subject of gas and liquid flow through solid porous 
bodies has been extensively studied. Flow of this type 
occurs in units designed for diffusional separations. Experi- 
mentalists interested in collecting design information for these 
units have presented a familiar correlation for countercurrent 
flow of gas and liquid through packings in which the pressure 
gradient in the gas phase is plotted as a function of gas mass 
velocity with liquid mass velocity as a parameter“). A number 
of approaches to single phase flow through porous bodies have 
been applied to two phase flow’). 

This work is a continuation of an earlier investigation 
into the nature of co-current and countercurient flow of gas and 
liquid through glass fibre beds. The effects of variables such as 
bulk density of the glass fibres, and the type of fibre are also 
inv estigated. The extremely fine nature of the bodies forming 
the solid bed strongly influences the hydraulic behavior and 
the results are significantly different from other reported studies, 


even though one or two have been concerned with the flow of 


two phases through glass fibres“. The conclusion from previous 
work is that the pressure drop across a given bed is defined by 
the liquid and gas flow rates. The work reported here shows 
that for some porous bodies at least, the pressure drop may be 
different for the same gas and liquid rates. 
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Apparatus and Procedure 


The arrangement of the apparatus is shown in the flow 
diagram Figure 1. The glass fibre bed was held in a 2-in. 1.D. 
glass column. Clean compressed air was delivered through a 
rotameter system to the top level of the fibre bed. The air 
pressure was controlled by a pressure reducing valve and a fine 
needle valve, so that the pressure of the gas entering the rotameter 
was close to atmospheric. Distilled water was drawn from a 
thermostatically controlled circulation system, measured by a 
series of rotameters and delivered to the top of the fibre bed 
through a liquid distributor. Pressure tappings at either end of 
the bed were connected to a differential vertical mneneeanet: 
The guaranteed accuracy of the rotameters was +2% of the 
indicated flow rates; the manometer was fitted with a vernier 
scale which allowed pressure differences as large as thirty inches 
of fluid to be read to 1/100-in. 

One of the most important operations was the preparation 
of the fibre bed. The arrangement of the glass fibres was 
standardized; each fibre lay within a plane normal to the main 
axis of the fibre bed and the fibres were randomly orientated 
within this plane. Discs of fibre matting of diameter equal to 
the internal diameter of the column were cut from a principal 
fibre mat, adjusted to a pre-determined weight and loaded into 
the glass column with each disc normal to the main axis of the 
column. The finished length of the fibre packing was constant 
for most of the experiments. The bulk density of the fibres 
was controlled for each bed. This was adjusted by increasing 
or decreasing the weight of the fibre mats in the standard length 
of column. Thus a ‘change in bulk density was achieved by 
increasing or decreasing the fibre density along the main axis 
of the column. The fibre density normal to the main column 
axis was not adjusted and remained constant. Variations in the 
bulk density of the fibres along the main axis of the column 
could arise if the packing was not carefully done and this distorted 
the hydraulic results. 

To assist in the definition of bed preparation, the variation 
of bed length with applied mechanical load was studied. It was 
found that the volume — load relation for the bed after proper 
consolidation was given by an expression of the form 

1 Vi — Vz, 
Pi = IS wo eee ee ee 

K Ve — V, 
where P; and P; are applied loads and V; and V2 are the bed 
volumes at conditions | and 2, K is a compressibility coefficient 
and V, is the volume of the bed as P+ (an extrapolated value). 
It was essential that the method of bed preparation gave a 
(P, V) value which followed this relation for this indicated a 
constant bulk density of fibre in the bed. If the mechanical load 
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Figure 1—Arrangement of apparatus. 


applied to the fibre bed was increased without attempting bed 
consolidation, unbalanced frictional forces acting at the ends of 
the fibre and the column wall produced a variation in the fibre 
bulk density which destroyed reproducibility of the hydraulic 
results. In this case the (P, V) values did not follow relation 1, 
but a marked hysteresis effect was found as the load was in- 
creased and then decreased. An even bulk density was assured 
if the bed was consolidated by alternately increasing and reducing 
the applied mechanical !oad, so that the (P, V) values were 

iven by Equation (1). At other (P, V) conditions the value 
of dV/0P depended upon the direction of the pressure change. 
This is a result of unbalanced frictional forces within the column. 

The bed when prepared was held in the column between two 
perforated plastic discs. The discs were supported between the 
pressure taps. The equipment was designed so that the fibre 
bed length could be adjusted and the whole column could be 
easily detached from its supports. The success of the method 
of bed preparation was accepted when the two phase pressure 
drops were reproducible for different beds of the same fibre 
packed in the same manner. 

The degree of reproducibility depended upon the nature of 
the fibre mat. If the mat showed little fibre aggregation, a 
representative sample of the fibre was easily obtained and the 
hydraulic results were reproducible within +2%. If fibre 
aggregation was prominent a small number of fibre mats did 
Not constitute a representative sample of the bulk fibre. Conse- 
quently the reproducibility limits were somewhat wider but in 
all cases better than +5%. 

Before making any pressure drop determinations, the fibre 
bed was treated with an aqueous solution of detergent to remove 
adsorbed impurities and to make sure that the fibre was properly 
wetted by the liquid. The detergent was removed by flushing 
the bed with water. 

The pressure drop across the bed at any set of flow conditions 
has been found to be dependent upon the previous flow history of 
the bed“. Consequently the starting condition was standardized 
and two principal starting conditions were used; liquid is allowed 
to flow through the bed for 15 minutes before the gas flow was 
directed through the bed (starting procedure 1) and the gas is 
allowed to flow through the water saturated bed for 15 minutes 
before liquid flow is directed through the bed (starting procedure 
2). Most of the determinations were made with starting pro- 
cedure 1. 

The liquid hold up in the bed was measured by removing 
the glass column from its supports and weighing. Because the 
liquid and gas flows are co-current, the liquid weight measured 
was the dynamic hold up. There was no drainage of liquid 
from the bed when the liquid and gas flows were shut off. 


Results 


The experimental results are shown in the series of graphs 
Figures 2-15. The numerical results are given more fully in 
Reference 5. Figure 2 shows the change of pressure difference 
across the bed with gas flow rate at a fixed liquid rate. The 
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Figure 2—The effect of gas flow rate upon pressure drop 
at fixed liquid rate. 


line ABC is obtained when gas flow is increased from starting 
condition 1 i.e., the gas flow is started when the liquid has been 
flowing for some time. The points A, B and C mark the maximum 
gas velocities in these experiments. From these points the gas 
velocity was reduced and the pressure difference — gas flow 
relation follows the lines AX, BY and CZ. Thus the operation 
of increasing the gas flow is not reversible. The graph also 
shows that the lowest pressure differences at fixed gas and 
liquid rates, are obtained if there has previously been a very 
high rate of gas flow. The lowest line on Figure 2 is the pressure 
difference — gas flow relation for the dry bed. 

Curves for this type were reproducible for the same bed 
and for other identically prepared beds. Thus there was no 
rearrangement of the fibres and, in fact, the fibre packing 
density was so high, that the fluctuations in gas pressure during 
an experiment did not affect the packing density. 

Points D, E and F record the pressure differences at the 
corresponding gas flows for starting procedure 2. The relevant 
liquid hold up measurements are shown in Figure 3. This graph 
shows that liquid hold up falls off as pas flow is increased, but 
Figure 5 shows that the irreversible fall in pressure difference 
is not matched by a comparable change in liquid hold up. In 
fact the change in liquid hold up was reversible with gas flow. 
The hold up values were independent of the previous gas flow 
rates, and the hold up — gas flow rate relation was single valued. 

Figure 4 shows the change in pressure differences at different 
gas rates for starting condition 1 with liquid rate shown as a 
parameter. Both liquid hold up and pressure difference increase 
with liquid rate. 

The effect of a change in bulk density is shown in Figures 6 
and 7 for starting condition 1. Both pressure difference and hold 
up are increased by an increase in fibre bulk density. 

‘ Figure 8 shows the experimental results for three types of 
bre. 


Table 1 shows the character of each type of fibre under test. 


The dispersion of the fibre bed refers to the manner in 
which the fibres are arranged relative to each other in a plane 
normal to the main axis of the bed. A well dispersed bed shows 
fibres which are randomly arranged without aggregation. Figure 
9 shows a photograph of a well dispersed bed. This may be 
compared with Figure 10, which clearly indicates fibre aggrega- 
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11-391 7.43 X 10-4 in 11.0 X 10-5 in. Poor 
11-390 8.53 X 10-4 in 9.0 X 10-5 in. Poor 
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Figure 4—The effect of liquid rate upon pressure drop. 





rr ~T T I 
FBRE NO, 3950 
BED LENGTH -084 inches 
BULK OENSITY A 


24 A+ 9°47 Ibs/tt® 
B= 8°60 
22 C= 671 







e | 
LIQUID FLOWRATE =529 Ibs/hy, ft 





$——$ ++ +—__ + 

| 

22 4 
“a 

4 

4 





LIQUID FLOWRATE =176 Ibs/hr ft® 








L 1 ce 1 L 1 
0 20 40 60 80 100 120 


GAS FLOWRATE Ibs/hr. ft® 


Figure 6—The effect of packing density upon pressure 
drop. 


tion and preferred orientation. Such a bed is poorly dispersed. 
Evidently the diameter and dispersion of the fibre greatly affect 
the hy draulic characteristics. All of the fibres, however, showed 
the type of behavior illustrated in Figure 2. The general 
effects of bulk density and liquid rate shown in Fi igures 4 and 6 
were confirmed for all types. 


Discussion 


The pore structure of the dry fibre bed is a continuous, 
twisting complex of interconnected channels bounded by the 
cylindrical walls of the individual fibres. Because of the method 
of bed preparation, the size distribution and arrangement of 
fibres in a small (but not differential) volume of the ‘bed is the 
same as any other small volume provided the bed dispersion is 
good. In two phase flow of air and water, the fibre is prefer- 
entially wetted by the liquid and the flow path of the liquid 
embraces the fibre and forms the continuous phase. The gas 
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Figure 3—The effect of a change in starting procedure 
upon liquid holdup. 
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Figure 5—The effect of liquid rate upon liquid hold up. 


flows through segregated channels in the bed and is the dispersed 

phase, bounded by liquid and solid fibres. The channels. bounded 

by solid in the bed are small, hence capillary forces are important 

at an interface, and as the fibre is wetted by the liquid, there is 

an excess of pressure in the gas phase over the liquid phase. 
The governing relation is 


p=T re ae Se Teer (2) 


where p is the pressure excess, T is the surface tension, and 7 
and rz are the principal radii of curvature at the elementary 
section of the interface considered. In this study, the gas phase 
pressure only was measured and we have no direct knowledge 
of static pressures in the liquid phase because of the capillary 
effect. The balance between the capillary forces and the fluid 
pressures determines the flow volume occupied by the gas and 
liquid phases. If the pore size is small, gas may be locked out 
of sections of the bed and can only advance into these sections 
if the gas phase static pressure is increased. 

For a small length of bed we can assume that the static 
pressure difference between the two phases is constant. Also 
as we have arranged geometrical similarity along the length of 
the bed, the distributions of flow for each phase are similar and 


211 










LIQUID FLOWRATE = 176 Ibs/ he, f1> 
BULK DENSITY 

A= 9°47 bs/tf 

B= 8-60 


Ce 67! 


WATER 


FRACTION OF BED VOLUME OCCUPIED BY 








° 20 40 60 80 100, 120 
GAS FLOWRATE  Ibs/hr, ft? 


Figure 7—The effect of packing density upon liquid hold 
up. 


random. The total volumetric flow rates are the same for each 
plane normal to the main axis of the bed. This suggests that 
the hydraulic characteristics may be studied through a model of 
continuous, irregular, interconnected channels whose walls are 
parallel to the main axis of the bed. The walls of the channel 
are formed from the supporting structures of the fibres and the 
interconnections are filled with liquid. The position of the 
liquid — gas interfaces are maintained by a balance between the 
capillary and pressure forces. 

Initially the fibres are wetted and a number of starting 
conditions are possible. We have defined two starting conditions 
and condition 1 probably corresponds to the bed voids filled with 
liquid. For the bed as a whole, gas flow can only be established 
in flow channels if the pressure difference across the bed exceeds 
the surface tension forces in the channel. 

If we use our model and apply it to the whole bed, the condi- 
tion for this is: 


R 
AP.A > T.Ror AP > ra: ease bh (3) 


where AP is the gas phase pressure difference, T is the surface 
tension, A is the cross-sectional area of the channel and R is 
the perimeter of the channel. The quantity A/R has the 
dimensions of a length and may be thought of as an ‘effective’ 
or ‘hydraulic’ radius. From this relation it can be seen that 


the larger channels will admit gas flow at lower pressure 
differences than the smaller channels. As the pressure difference 





Figure 9—Structure of well dispersed 
fibre bed. 
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Figure 10 — Structure of poorly dis- 
persed fibre bed. 
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Figure 8—The effect of fibre type upon pressure drop. 


across the bed is increased, a greater number of channels will 
admit gas flow and consequently the resistance of the bed to 
gas flow will decrease. 


We have assumed that the pressure excess of the gas over 
the liquid phase is constant, independent of length, for a small 
length of bed. Suppose that the constant value is increased — for 
example by increasing the static pressure of the gas phase and 
the flow of gas. The position of the gas — liquid interface at 
each point is adjusted until a new position of balance between 
the capillary and pressure forces is reached. Figure 11 illustrates 
a possible situation. 


The line AB marks the initial position of the gas — solid 
interface before the gas flow is increased. As the gas phase 
pressure grows larger, the interface position is adjusted so that 
the curvature (1/r; + 1/72) is increased. This means that the 
interface retreats to a position where the interfacial area is 
smaller. From the diagram it can be seen that this is not accom- 
plished until a portion of wider section has been traversed and 
a balance is reached at the new position, CD. If the gas phase 
static pressure now is reduced to its former value, the interface 
will move back to EF and not to AB. Clearly to return to 
position AB the gas pressure would have to be lowered to a 
value sufficient to move the interface past the widest part of 
the flow section. Thus the increase in gas flow rate has led to 
an irreversible change in the flow distribution of the two phases. 

The Navier-Stokes equations describe the flow of gas through 
the bed and for the simplified case of the model the gas flow may 
be assumed to take place in one direction and two dimensions. 
When the flow is constant at steady state these equations 
reduce to: 

Ou O% 1 0P 


: mom OL. 4 
xt * dy aie 0 (4) 





Figure 11—An irreversible change in 
gas and liquid flow distribution. 
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Figure 14—The effect of bed length upon pressure drop. 


where w is the gas velocity, u is the gas viscosity, P is the gas 
phase pressure and x, y and z are space co-ordinates. 

Throughout the experimental work, the rate of gas volu- 
metric flow was much greater than the rate of liquid volumetric 
flow, and as the liquid hold up determination showed that the 
volume of the bed occupied by liquid was greater than that 
occupied by gas, it follows that the gas phase velocities in the 
bed were generally much higher than the liquid phase velocities. 
The boundary conditions for each flow channel are vu, = u at 
the liquid gas interfaces and u, = 0 at the solid boundary. As 
the liquid phase velocities are much lower than the gas, we can 
assume that a general boundary condition is given by u, = Oat 
the solid or liquid boundary of the gas phase. 


The differential equation is a form of Poisson’s equation if 


dP/dz is constant. With the condition of zero velocity at the 
boundary, the solution to this equation which gives u as a 
function of x and y may be studied by analogy with Prandrtl’s 
loaded membrane. Thus the displacement of a membrane 
secured at the boundary and subjected to a pressure difference 
will give the form of the relation. 

The solution gives “, = 0 at the boundary and a maximum 
at points most distant from the boundary. 
velocity is a maximum at the boundary and falls to zero at the 
points of velocity maxima. The general nature of the solution 
may be used to ‘illustrate an important relation concerning the 
resistance to flow of the model. 

An arbitrary gas flow section is shown in Figure 12. The 
shear stress 7 at the boundary of the flow section is 


T = wgradu... - (5) 


where wu is the fluid viscosity and grad w is the derivative of 


velocity measured normal to a line of equal velocity on the 

flow surface. 

increasing velocity. If we equate pressure and viscous forces 

over the flow channel during gas flow 

‘A B 
megrad udR + H 


B 4 


AP.A =H mw grad udR (6) 
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Figure 15—The effect of bed length upon mean pressure 
gradient. 


where A is the cross-sectional area of the flow channel, AP is 
the pressure difference for a length of channel H and R is the 
B 


perimeter of the boundary. H | ym grad wdR is the shear force 
A 

acting at the short portion of the boundary marked AB. Suppose 
that the gas phase static pressure is increased and reduced to its 
former value in the manner we have described. As a result of 
this the interface position AB is now moved to CD. If we again 
equate pressure and viscous forces over the original channel: 

A B Ou 


uegradudR+H} up dR. (7) 
B a \o% 


where “(0u/d@) is the shear stress normal to the boundary 
AB. The left hand side of Equations (6) and (7) remains 
unchanged. The last term on the right hand side is the shear 
force acting over the section AB. This is reduced when the 
boundary is extended since the velocity along AB is increased 
from zero and hence du/d@ along the boundary is reduced. 
Thus the shear force represented by ‘the last term in Equation (7) 
is smaller than the last term in Equation (6). The shear force 
represented by the first term in the RHS of Equation (7) is 
larger than the first term on the RHS of Equation (6). This 
shows that the velocities in the region of the solid section of 
the boundary have increased. The total gas flow in the original 
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section is increased at the same — gradient by extending 
the boundary in this arbitrary manner 

This may be expressed: If the ends of a viscous fluid 
flowing i in a channel of constant cross-sectional area is extended 
in any way, the resistance to flow of the original channel is 
decreased. For a single gas flow channel an increase in gas flow 
rate followed by a decrease to the same value will lower the 
resistance of the original channel to gas flow as some of the 
extension of the flow channel is retained. It follows that the 

as flow resistance for unit gas flow area in the fibre bed is 
reduced by an increase in the gas flow rate followed by a decrease. 

The extension of the gas flow channels reduces the liquid hold 
up. There is however, a counteracting effect; the reduction in 
pressure drop across the bed corresponding to the reduced flow 
resistance will close off some of the gas flow channels according 
to Equation (3). These channels, active before the increase in 
flow rate are filled with liquid and this will increase the liquid 
hold up. Under the experimental conditions of this investigation 
these effects were found to be counterbalancing and the liquid 
hold up — gas flow relation was reversible. 

From the study of the flow equation and the flow behavior 
in the bed, it follows that an increase from, followed by a 
decrease to the same gas flow rate will extend some of the larger 
gas flow channels. As a result the resistance of the fibre bed 
to gas flow will fall even though the porosity (liquid hold up) 
is unchanged. 

Figure 13 shows the change in flow structure of the bed as the 
gas rate is first increased and then decreased at a fixed liquid 
rate. As the gas rate is increased, more channels become open 
to gas flow in accordance with Equation (3) and existing flow 
channels are extended. As gas flow is reduced, the extended 
channels retain some of their extension, thus the resistance of 
the bed to gas flow is lowered and the pressure difference across 
the bed is reduced. Consequently, fewer channels are open to 
gas flow. The aggregation of gas flow channels is shown in the 
marked reduction of resistance to gas flow illustrated in Figure 2. 


In the development we assumed that the static pressure differ- 
ence between the gas and liquid phases is constant for a short 
length of bed. It is of interest to see if this is true for the total 
length of bed. Figure 14 shows the pressure difference across 
the bed plotted as a function of bed length. Figure 15 shows the 
mean pressure gradient plotted as a function of bed length. The 
gas phase mean pressure gradient falls off at equal gas and liquid 
rates as the bed length is increased. Thus the fall in gas phase 
static pressure across the bed reduces the static pressure differ- 
ence between the two phases. Consequently the liquid hold up 
is larger at the lower levels of the column and as a result the 
relation between gas phase static pressure and bed length is non 
linear. The effect is more marked at lower gas rates. Figures 3 
and 5 show that the change in liquid hold up is most pronounced 
at the lower gas rates. The effects are evidently the same. 


Irreversibility of the Change in Flow Structure 

The study of this flow system has shown that an increase 
in the gas flow rate at a fixed liquid rate, causes an irreversible 
change in the flow structure of the bed. Consider a cycle of 
flow changes: the gas flow rate is increased at a fixed liquid 


g = gas 
l 


flow rate, the liquid rate is then decreased to a lower value 
maintaining the gas flow rate constant, the gas flow rate is 
reduced to its first value, and finally the liquid rate is increased 
until the original pressure difference across the bed is attained. 
If the flow changes do not cause an irreversible change in the 
flow structure of the bed, the liquid rate at this point is the same 
as the initial value: 
,, dl = 


where @,, indicates that the integration is carried round a 
complete “cycle of pressure difference- “gas flow conditions. If, 
however, there is an irreversible change in the flow structure of 
the bed: 


Combining these two expressions, for any cycle of changes on 
the P — G diagram 


where the equality sign oie to reversible changes, and the 
inequality sign applies to irreversible changes. This integral is 
thus a measure of the irreversibility of the operation. 
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Nomenclature 
A = area of cross-section 
G gas phase mass velocity 
H length 
K = compressibility coefficient 
i liquid phase mass velocity 
ln natural logarithm 

static pressure excess of gas phase over liquid 
ar load, pressure 
R perimeter 
radius of curvature 
sutface tension 
velocity 
volume 
space co-ordinates 


<2 a 
. 3 


” ” 


9» ” 
difference quantity 
viscositv 

shear stress 


WE pr*e a 


Subscripts 
liquid 
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Pelletization and Diffusibility Studies 
on Silver Catalysts’ 


C. H. AMBERG* 


Powders of silver catalyst of low calcium content 
and an 8.5% calcium in silver alloy were pelletized 
at different compacting stresses. The diffusibility of 
each pellet was measured, using hydrogen, and cor- 
related with internal void volume. The surface area 
of the powders and of some representative pellets was 
obtained by nitrogen adsorption. All these parameters 
were used in an attempt to interpret the changes 
taking place with decreasing porosity. With the help 
of two plausible assumptions it was concluded that 
the tortuosity remained constant over the range 
examined. The figures obtained have limited useful- 
ness for the prediction of approximate pellet diffu- 
sibilities from a knowledge of the compacting stresses 
applied to similar metal powders. 


9 heterogeneous catalysis a knowledge of the internal structure 
of a catalyst is frequently desirable in order to evaluate the 
contribution of diffusion to the rates of transport of different 
chemical species to and from the active surface”. In this 
laboratory porous sprayed silver films have been employed in 
the catalytic oxidation of ethylene. In the hope of approximating 
their structure in a controlled fashion, several series of pellets 
of varying porosity were prepared from a powdered Ag-Ca 
alloy and from silver catalysts originating from this and similar 
alloys. The diffusibility of each pellet was measured using 
hydrogen; this was correlated with the internal void volume in 
an attempt to gain some insight i into the highly complex geometry 
of the compact. In some instances the available surface : area was 
determined by nitrogen adsorption. We have also compared 
compacting pressures with the void volumes so that one might 
be able to prepare pellets of predictable diffusibility. However, 
as will be seen this relationship lacks generality; only very 
approximate diffusibilities can be predicted unless extensive 
experimental information about a given catalyst powder is 
available to begin with. 


Experimental 
The apparatus used was essentially the same as that described 


by Weisz and Prater which, in turn, is a simplified version of 


an earlier one devised by Wicke and Kallenbach®. It contains 
a sample pellet exposed on one side to a stream of hydrogen, on 
the other to one of nitrogen. Flow rates were such that there 
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was zero hydrostatic pressure difference across the pellet, while 
maintaining a partial pressure drop of one atmosphere. 

The concentration of hydrogen diffusing through the sample 
into the nitrogen stream was measured by a ‘thermal conductivity 
cell (Gow-Mac Co., Model 9285) immersed in a bath thermo- 
regulated at 30.0°C. Flow rates, including those used for 
calibrating the conductivity cell, were read from differential 
manometers. Pure nitrogen was employed as the reference gas 
in the thermal conductivity cell before being passed into the 
diffusion cell. 

Pellets were prepared in a simple right cylindrical punch- 
and-die assembly. Pressure was applied from opposite directions 
to the two punches (diameter approximately 1 cm.) for one-half 
minute, no attempt being made to study time effects with the 
continued application of pressure. The pellets, usually 2-4 mm. 
thick, were mounted in Tygon tubing as described by Weisz 
and Prater“. 

Diffusion measurements were carried out at room tempera- 
ture and corrected to 25°C. Each figure quoted comprises the 
average of about five determinations over a wide range of 
nitrogen flow rates. When the product of nitrogen flow rate 
and hy drogen concentration remained constant within 1%, this 
quantity was used to calculate the effective diffusion coefficient 
Dg and hence the diffusibility 6 = D./D. Failure to attain 
constant values was taken to mean that the pellet was defective 
or imperfectly mounted. 

The main materials employed were an Ag-Ca alloy contain- 
ing 8.5% Ca, a silver catalyst prepared from this alloy by steam 
and acid treatment“ and a similarly prepared catalyst. The 
latter is designated as No. 1 (Figures 1 and 2), while the alloy- 
derived cataly st and its parent alloy are Nos. 2 and 3, respec- 
tively. The two catalysts differed somewhat in their residual 
Ca content, which amounted to 0.10% Ca for No. 1 and 0.06% 
Ca for No. 2. Spectrographic analyses of trace elements were 
similar, except that“No. 1 showed slightly higher Al and slightly 
lower Fe content than did No. 2. B.E.T. nitrogen surface areas 
were 0.4; m.?/g. for No. | and 0.2 22 m.*/g. for No. 2, while 
the alloy gave a value of 1.50 m.?/g. A variety of mesh sizes 
(Tyler) and size ranges were used; these will be found indicated 
in the Fi igures. 

Additional materials used in a few separate experiments will 
be described along with the pertinent results. 


Results 

In Figure 1 are shown the correlations between the axial 
stress exerted on the pellets during compacting, plot.ed loga- 
rithmically, and the void fraction. It is customary to report 
compacting stresses as pressures p = F/A; we have followed 
Duffield and Grootenhuis® in employing what is probably a 
more realistic quantity, namely the axial stress s = FTp/W. 
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Figure 1—Void fraction vs. axial stress. 1 and 2, silver 

catalysts; 3, Ag-Ca alloy. Mesh size ranges: 35-80, O; 48- 

100, ™ ; 80-100, (J; 100-150, /\; 100-200, © ; 150-200, 
@; 200-250, MJ. Insert: see text. 


In contrast to p, which does not differentiate between solid and 
void, s represents an average stress on the solid matter only. 
In calculating s, the “dry friction” which causes a decrease in 
axial stress fond thus pellet density) as one proceeds axially 
from the faces to the interior of a pellet was neglected. The 
error thus introduced was probably not larger than 3%, 
Possibly similar errors were introduced by changes i in packing 
density in a radial direction. Thus a photomicrograph taken 
near the rim of a lightly compressed silver pellet, sectioned 
parallel to its faces, ‘showed a 2% area of somewhat higher 
porosity. This is never evident from the outside appearance of 
a pellet as the rims look always particularly smooth. The void 
fraction @ was calculated from the external pellet dimensions 
and the density of the solid. It thus includes both open and 
totally enclosed, i.e., inaccessible, void spaces. 


All three types of metal powders gave different compacting 
curves, which were linear over the regions examined. Curves 
1 and 2 appear to be separated by an approximately constant 
void fraction increment of 0.1. The alloy powder, curve 3, 
was less easily compressible; moreover it was not possible to 
form coherent compacts with void fractions in excess of 6 = 0.3. 
As indicated in the figures, each type of powder was used over 
a fairly wide range of particle sizes. Void space formation for 
the two catalysts appeared to show no evidence of dependence 
on particle size or, within the very limited number of tests, 
extent of size range. In the case of the alloy it was particularly 
difficult to form good pellets from coarse particles, which may 
have given rise to errors as reflected by the spread in points. 
On the other hand it is possible that the effect was real. Since 
the alloy powders had been prepared by crushing and several 
stages of grinding, the largest particles had experienced the 
least amount of shear. Thus different size ranges may well 
have exhibited variations in surface structure which, in turn, 
became less prominent after the leaching process. 

The framed area about curve 2 is redrawn in the insert to 
Figure 1 with the ordinate enlarged 21% times. All points save 
the lowest were obtained with catalysts similar to No. 2 and 
serve to illustrate the possible effect of composition and pre- 
treatment on pellet formation at a constant compacting pressure 
of 2830 kg./cm.* Points a and b are duplicates for a catalyst 
somewhat higher in residual Ca(0.33%) and presumed to have 
a surface area higher ~~. No. 2. Point ¢ represents a catalyst 
prepared from alloy No. 3, but by a method resulting in more 
than double the surface area (0.48 m.? /g.) at the same residual 
Ca content (0.06%). Catalyst No. 2 was rechecked (d) and 
subjected to heat treatments (¢). Identical values were obtained 
after 6 minutes at 425° or 10 hours at 350°C, This is of interest 


216 


DIFFUSIBILITY 8 








Na i 
00 Ol “sas aa 05 0607 00 01 0203 
VOID FRACTION @ 
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Figure 3—B.E.T. surface area relative to that of loose 
powder vs. accessible void fraction, for catalyst No. 2 and 
alloy. 


in view of the equal activities in ethylene oxidation achieved in 
this laboratory after comparable pretreatments of a similar 
catalyst. Finally, point f was determined with a pure commer- 
cial silver powder. 

For purposes of orientation commercial powders of Cu, Ni, 
ZnO and A1,O3 were pelletized and the results indicated in the 
diagrams. The Cu, particularly, consisted of smooth, nearly 
spherical particles and the Ni particles, while not quite so 
regular in shape, still presented a far smoother appearance under 
the microscope than the silver catalyst. 

The diffusibility values as a function of void fraction are 
shown in Figure 2, the straight lines having been computed by 
the method of least squares. Curves No. 1 and 2 were similar 
enough so that they were treated as arising from one and the 
same catalyst. In the case of the alloy the three values for the 35- 
80 mesh powder were omitted in the calculation. In both cases the 
extrapolated lines intersect the abscissa at finite values of @ = @, 
which are, respectively, 0.095 and 0.040 for the catalysts and 
the alloy. 

Following the method of Weisz and Prater we have 
conducted a number of forced flow experiments with air and 
have found no evidence for Knudsen flow. This is to be expected 
also on the basis of pore sizes (estimated from our photomicro- 
graphs as well as from initial particle diameters), which lie well 
within the region of Poiseuille flow at atmospheric pressure”. 
In this connection it should be mentioned that the extrapolations 
to zero diffusibility are likely to be fictitious, although for the 
catalysts departures from the straight line occur only within a 
small region. A forced flow experiment with the pellet charac- 
terized by the lowest point was unfortunately inconclusive 
because of the extremely low rate of flow involved. 

B.E.T. nitrogen surface areas, 2, were determined for catalyst 
No. 2 and for the alloy at several stages of compression. In 
Figure 3 these are plotted as fractions of the value of the 
uncompacted powder, 2), against @ — @,. The void fraction of 
the powder bed was obtained from the dimensions of the bed 
after tapping to constant height. 
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Discussion 


The relationship between pelletizing pressure and total void 
space of the compact formed has been extensively examined by 
Duffield and Grootenhuis for spherical copper particles“, where 
equations proposed by other workers were also critically dis- 
cussed. In view of their findings no general significance should 
be attached to the approximate linearity exhibited by our curves 
in Figure 1. However, for the purpose of preparing pellets of a 
given void volume and diffusibility it should prove useful in an 
empirical fashion once initial data have been determined for a 
given powder. 


Duffield and Grootenhuis noted that the behavior under 
compression was independent of particle size except through 
its influence on the initial packing which, in turn, governed the 
absolute position of their pressure-volume curves. Within the 
limits of accuracy of the present investigation there was again 
no dependence on particle size. Moreover, it may be presumed 
that initial packing, in this case strongly determined by particle 
shapes and, at least for the catalysts, by easy agglomeration 
under low stress, again played a decisive part in fixing the 
absolute value of the void space at a given axial stress. Indeed, 
the greater roughness of catalyst No. 1 1, also observed visually 
when compared to No. 2 under a low-power microscope would, 
if this picture is correct, provide the basis for the observed 
higher void volume when the particles become “locked” in 
position upon the initial application of pressure. It is unlikely 
that the greater volume can have been caused by significant 
differences in sack pores or totally enclosed voids, since these 
would have given rise to differences in diffusibility between 
Nos. 1 and 2 at a given total void fraction. This was clearly 
not observed (Figure 2). 


The inset to Figure 1 shows the influence of Ca content and 
surface structure on void formation for a constant compacting 
pressure, with particle geometry being the predominant factor 
at least in the range 0 — 0.3% Ca. It might be added that owing 
to the method of computing @, reproducibility is only about 
+3% in these cases. The total pellet volume per gram on the 
other hand was reproducible to within +1% and might well 
be used to good advantage for a rapid characterization of the 
macrostructure of catalyst powders. 


At the higher Ca content of the alloy (curve 3) conditions 
were markedly different. The material is more brittle and, 
although we have worked with freshly prepared samples, there 
is always the possibility of CaO formation. This would provide 
a barrier to surface diffusion and the growth of inter-particle 
necks, all necessary elements of the sintering process (see, ¢.g., 
®), Both factors would be reflected in the smaller slope 
of curve 3. The spread in points with particle size has already 
been referred to; the curve as drawn is biassed towards the 
smaller size ranges. 


It is interesting to note that the microscopically fairly smooth 
particles of pure Ni, Cu and Ag compressed at 5230 kg./cm.? 
(Ag not shown in Figure 1) form pellets with void fractions 
correlating approximately linearly with 7}, the melting point 
of the pure metals in °K., and very roughly in line with their 
respective shear moduli (see Table 1). The value of 0.5 77, is 
often taken to signify the approximate temperature at which 
two like metal surfaces form adhesive bonds under normal load. 
Rowe has measured the 7, factors to be 0.45 and 0.40, 
respectively, for clean Ni and ‘Ag cylinders in a high vacuum. 
Our correlation could be interpreted in terms of the significance 
of surface mobility and malleability in the processes under dis- 
cussion, but data on further metals are obviously required to 
taise this out of the realm of speculation. 

The values for ZnO and Al,O; shown in Figure | may be 
explained on the basis of much higher resistance to sintering 
and deformation, as well as microporosity in the case of ALOs. 
Thus, Taylor and Amberg found relatively little loss in surface 
area”) on pressing high-area ZnO powders under loads up to 
ten times those employed in the present investigation", 
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TABLE 1 


COMPARISON OF VOID FRACTION WITH MELTING POINT AND SHEAR 
MODULUS. 3. COMPACTING PRESSURE: 5230 kg/cm?. 


| 
| 





is aa Shear Mod. 

c 0.5 Tn°K dyne/cm.? 0 
Ag | 617 | 2.9 x 10" 0.112 
Cu | 678 4.6 X 10" 0.130 
Ni | 864 7.7 X 10" 0.256 


The meaning of the diffusibility 6 in terms of idealized pore 
structures has been discussed fully by Zimens“* and Wheeler“. 
Discrepancies between and @ are customarily ascribed to 
tortuosity, expressed as x = 6/8. Hoogschagen“” has quoted 
values of x for various ocihicieal catalyst powders, which were 
of a magnitude similar to ours (Figure 2). Apart from the contri- 
bution of the true tortuosity of the path of the diffusing fluid, x 
when computed in the above manner must also include contribu- 
tions from void spaces not participating in fluid transfer. These 
could be either totally enclosed voids or “‘stagnant’”’ volumes 
across which the partial pressure gradient of the diffusing species 
is zero. 

If we accept the evidence for Poiseuille flow in our compacts, 
it would follow from the work of Wheeler and others that 
the approximate linear relationship exhibited for the curves 
shown in Figure 2 indicates an essentially uniform decrease in 
pore sizes without change in pore number or geometry over 
the linear region. This assumes that there is no fortuitous 
occurrence of counteracting effects (e.g., an increase in the 
size of large pores at the expense of smaller ones and decrease 
in pore number) which is considered unlikely because of the 
similar behavior of the physically different systems 1 + 2 and 3. 
It is also compatible with the measurements of the internal 
surface area discussed below. We shall now assume as a first 
approximation that the “stagnant’’ void volume is negligible 
and that the blocked void fraction @, is a constant throughout 
the linear region. It is difficult to assess the validity of these 
assumptions. No experimental method is available to us for 
measuring the first quantity and helium density determinations 
have failed, for lack of the required precision with the available 
samples, to define the second. For the blocked void fraction to 
be constant with decreasing pellet volume per unit weight, the 
blocked void volume must decrease, which seems reasonable, 
especially since one may expect the immediate creation of 
blocked spaces on first forming a coherent compact of “open” 
structure. If this were not so, one should have expected marked 
curvature in the diffusibility curve, or even a kink at the onset 
of blocked pore formation. However, approximate linear be- 
havior of @, with 6 other than constancy cannot be discounted. 

In the light of the foregoing arbitrary assumptions we may 
express the graphs in Figure 2 by the relation 


6 = x(0 — 8,), 


where x is now a representation of the actual tortuosity. In 
previous studies“ it had always been caluclated from separate 
pairs of values, which in the present case would have resulted 
in a tortuosity increase (decreasing value of x) with increasing 
compression; the above expression indicates maintenance of 
essentially constant tortuosity over the whole range. The 
least square values of x so obtained were 0.17 + 0.013 and 
0.20 + 0.01, for the catalyst and the alloy respectively. 

We have already commented on the relative lack of depend- 
ence of the pressure-volume curves on particle size for a given 
type of powder. Figure 2 lends support to the idea that particles 
initially interlock forming agglomerates which are of similar 
geometry for different particle sizes because individual particles 
have lost their identity. Indeed, the similarity of curves 1 and 
2 would suggest that the geometry of the two catalysts is essen- 
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tially the same, even though, for reasons mentioned earlier, 
different loads had to be applied to create pellets of the same 
total void fraction. 

With the additional assistance of Figure 3 we shall attempt 
to describe the processes involved in pellet formation. Perhaps 
the simplest case to interpret is that of the alloy. From measure- 
ments of the available surface area it appears that there was a 
tenfold decrease in area on first compressing the powder into a 
coherent compact of the maximum possible void volume. The 
comparatively small value of 8, can only mean that this decrease 
was caused either by removal of roughness through shear forces) 
and the formation of large interparticle contact areas (i.c., flatten- 
ing of particles), or possibly the closure of internal fissures, rather 
than by entrapping large voids. One might expect the shearing 
process to be predominant at this stage, since one could expect 
particle flattening to be accompanied by substantial densification 
of the pellet, w hereas there was in fact little change in void 
fraction between the loose powder and the first compact. On 
further compression the available surface area remained constant 
while the void fraction decreased. At the same time the diffus- 
ibility decreased linearly with void fraction. These facts rule 
out any effective destruction of pores, either by closing their 
outlets or by complete densification, and again point to mainten- 
ance of the original geometry. One may visualize this, for 
instance, as a narrowing down of flat-walled pores. 


The situation is somewhat more complicated with the 
catalyst. For reasons referred to above it is understandable 
that coherent compacts could be formed at pressures five times 
as low as for the alloy, resulting in a total void fraction as high 
as 0.6. The surface area did not materially decrease on reducing 
the void fraction to ~ 0.4, after which it decreased rapidly. 
Thus the initially wide open structure was first condensed by 
increasing the number of inter-particle weld-necks or the 
contact area of the existing ones, without however changing 
the pore geometry to any ‘significant extent, as evidenced by 
the approximately constant surface area and the linear decrease 
in diffusibility. Further compression occurred with a concomit- 
ant decrease in surface area, probably involving surface or 
bulk diffusion processes and plastic deformation, but still did 
not involve a change in flow pattern. Possibly deo a surface 

“‘smoothing”’ by shear as described for the alloy may have been 
involved, although it should be noted that 2, for this catalyst 
was only 0.2 m.?/g., whereas it amounted to 1.5 m.?/g. for the 
alloy. This meant that the catalyst presented a much ‘smoother 
surface to begin with and probably contained far fewer mechani- 
cally weak protuberances than did the alloy. 
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Nomenclature 
A = area of pellet face, cm. 
D = diffusion coefficient of Hz in No, at 25°C., cm.?/sec. 


Des = effective diffusion coefficient in pellet, at 25°C., cm.2/sec. 
’ = applied load, kg. 

Pp = compacting pressure, kg./cm.? 

= axial stress, kg./cm.? 

= pellet thickness, cm. 

= melting point, °K. 

pellet weight, g. 

diffusibility, dimensionless 

= total void fraction of pellet, dimensionless 

= extrapolated void fraction of pellet at 6 = 0, dimensionless 

= density of solid, g./cm.® 

= B.E.T. nitrogen surface area of pellet, m.?/g. 

= B.E.T. nitrogen surface area of uncompacted powder, 

m.?2/g. 

tortuosity, 
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Chemical Engineers and Computers’ 


THEODORE J. WILLIAMS? 


Computers as tools for the chemical engineer have 
begun to revolutionize his concepts of reaction kinetic 
analysis, of plant design methods, and of plant control 
systems. The new field known as systems engineering 
has arisen to fulfill and exploit the promises which 
computers offer in these areas. This paper outlines 
the applications of analog and digital computers to 
chemical engineering problems, points out the specific 
advantages of each, and details the numbers and 
backgrounds of the operating personnel required. 
The newer field of application of computers to process 
control is reviewed in the light of progress to date 
and of present and future trends. 


an are now revolutionizing chemical engineering. As 
tools of the design engineer they can release him from the 
hours of drudgery formerly associated with the trial and error 
aspects of design. They permit him to try many, many more 
different and vastly more complex design concepts. As tools of 
the research engineer they can mechanize the data collection 
and data reduction aspects of his experiments. They permit 
him to try out very complex mathematical models of his system 


prior to experimentation so that only the minimum number of 


the most crucial experiments need be made. As tools of the 
plant production engineer they can remove the necessity for 
close monitoring of plant operation by human operators; they 
can smooth out plant operations; and they can assure that the 
plant is always operating at its economic production level. 


In fact, the impact of the computer on chemical engineering 
has been so great that a new branch of chemical engineering, 
known variously as systems engineering, or dynamic analysis, 
or applied mathematics, or by many other terms depending upon 
the group involved, has come into being. To implement the 
uses of computers: listed above, systems engineering has taken 
the charter given in Table 1, ‘As a result of this and allied 
work many important changes have already taken place in 
chemical engineering studies. Some of these are: 


(1) The most striking changes have come in chemical reaction 
kinetics. Where older methods required the development of 
overall orders of reaction and overall rate constants, the newer 
methods along with computers permit the specification of each 
of the reaction steps and their individual rate constants. While 
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the resulting systems of equations are very large and the deter- 
mination of the correspondingly large ‘numbers of reaction 
coefficients can be quite difficult, computing machines make 
such work now practicable. 


(2) Methods are now being devised which promise to enable 
the determination of the phy sical properties of gaseous and 
liquid systems from readily made measurements of only a few 
of the phy sical properties ‘of the components of the sy "stem, 
Such methods will be of immeasurable aid to plant design studies 
when finally perfected. 


(3) Current experimental research along with the application 
of servomechanism techniques is rapidly establishing the dynamic 
behavior of process equipment when under transient conditions 
as when subject to upsets or during start up. 


Analog and digital computers each have specific advantages 
and disadv antages in relation to a particular problem which must 
be kept in mind in applying them. Table 2 summarizes this 
information. 


The short but spectacular history of computers, their ex- 
ponential increase in numbers, and capacity, and the fact that 
they have only been applied to chemical and petroleum industry 
problems to a high degree within the last five or six years 
means that most of the major discoveries concerning their use 
in our industry are still ahead of us. 


TABLE 1 


THe Aims OF SYSTEMS ENGINEERING IN APPLYING 
COMPUTERS TO THE CHEMICAL ENGINEERING FIELD 


(1) In Plant Design 

(a) To optimize the process while in the design stage, so as to obtain the 
best possible combination of processing conditions and plant design to 
assure the maximum return on investment. 

(b) To use computer techniques in plant unit designs, especially where 
trial and error procedures are necessary, to obtain the best possible 
designs for the plant units involved. 

(c) To correlate process dynamics and process control studies directly 
with process unit design studies. 


(2) In Research and Process Development 

(a) To render as much service as possible to the research chemist in 
developing the nature and characteristics of the chemical reactions 
involved in a process. 

(b) To make the fullest possible use of computer simulation in the process 
development stage in order to reduce the experimentation necessary 
in proving the process and establishing the most desirable operating 
conditions and the data necessary for the design phase. 


(3 


In Plant Optimizations and Plant Control 


(a) To study the operation of each existing plant to determine what, if 
any, changes in their control systems or methods of operation would 
improve their return on investment. 

(b) Study each of the rewer techniques of plant control to determine their 
applicability to existing processes. 
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TABLE 2 


A COMPARISON OF THE APPLICATION OF ANALOG 
AND DiGITAL COMPUTERS TO CHEMICAL ENGINEERING PROBLEMS 


(1) The Advantages of Analogs in Comparison to Digital Computers. 


(2 


Analog Computer 
Simulates behavior of any system 
by action of easily manipulated and 
measured variables. 


Simulation is continuous, permits 
inclusion of concepts such as dis- 
tance, velocity, acceleration. 


Results presented as family of 
graphs of variation of dependent 
variable (same data obtained with 
recorder connected to process). 


Speed of problem solving is direct 
function of actual speed of physical 
system, is independent of size or 
complexity of system. 


Programmed so that parameter 
magnitudes are entered as settings 
on variable potentiometers; para- 
meter values can be changed at will. 


Digital Computer 


Performs arithmetic operations 
with numbers. 


Operates discontinuously; can only 
approximate higher order effects. 


Results presented as tables of 
numbers; must be plotted. 


Speed is direct function of problem 
complexity and size, relatively 
independent of operating speed of 
process. 


Flexibility not present without 
special programming precautions. 


) The Advantages of Digital Computers in Relation to Analog Computers 


Digital Computer 
Problem size reflected in computing 
time; no limit to size of problem it 
time is available. 


Useful in nearly all types of prob- 
Jems. Complex partial differential 
equations more amenable to solu- 
tion by digital equipment. 


Analog Computer 


Size of computer (number of com- 
puting components) determines size 
of problem which can be solved; 
simplifying assumptions often made 
to reduce complexity. 


Ssecifically designed for solution 
of ordinary differential equations; 
other problems — e.g., solution of 
simultaneous algebraic equations — 
made by trial and error. Less 
complex partial differential equa- 
tions solved if converted to ordinary 
differential equations and one vari- 
able assumed constant; family of 
solutions for different values of 
variable results. 


(3) The Relative Applicability of Each Type to Specific Problem Areas. 


| 
Application 


Analog 


Digital 


Equipment and 
Plant Design 


Instrumentation and 
Data Reduction 


Proc ess C ontrol 
Simulations | 
(Research) | 


Cc he mical Plant | 
Applications 


Research in Basic 
Processes 


Resez arch in Unit 


Plant Operations 
and Management 
Technique Studies 


Show great promise for 
working out relative 
sizes of plant equip- 
ment by solving 
dynamic equations of 
plant performance. 


Small, special purpose 
analogs important 
in converting and 
plotting single variable 
and directly correlated 
multiple variable data 
in form of graphs. 


Best tool in existence 
for studying resulting 
phenomenon; limita- 
tion on system com- 
plexity imposed by 
computer size. 

Definite promise as 
final control elements 
for units of highly 
automated chemical 


plants. 


Well suited to simula- 
tion of chemical reac- 
tion kinetics involved 
in chemical process 
development. 


Well suited for study 
of transient state such 
as heat transfer, distil- 
lation. 


Superior for details of 
mechanical design — 
e.g., stresses, sizing of 
members, costs — 
which are trial and 
error arithmetical 
problems. 


Flexibility and 
accuracy of digital 
computers better for 
control of whole plants. 


Better for determina- 
tion of kinetic para- 
meters by statistical 
data reduction where 
necessary. 

Better suited for 
steady-state studies 
which are usually 
arithmetical or statis- 
tical. 

Better here because 
statistics and arith- 
metical operations 
predominate. 


However, this in no way implies that the difficulties and 


research necessary to apply computers is at an end. 


As an 


example of the status of one branch of chemical engineering 
the field of reactor dynamics and reaction kinetics will be 


discussed. 
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Reactor Dynamic Analysis and Process Kinetics 


The most important class of chemical processing units from 
a dynamic behavior viewpoint and because of its ultimate 
importance to the chemical process involved is the chemical 
reactor. Therefore, the most important research task of the 
process control or systems engineer is to determine the mixing 
characteristics of these devices and the kinetic mechanisms of 
the chemical reactions taking place therein™. 


Mixing studies are important in order to define the degree 
of complexity to which analog or digital computer simulations 
of the reactor must be carried out. If a stirred tank reactor can 
be considered completely mixed only one stage needs to be 
simulated On the other hand, inadequate mixing may require 
the simulation of 2 or 4 or even more stages per reactor. Thus 
if a very complex kinetic model is required this may automatic- 
ally make the simulation impossible to ) wey out on available 
computing machines. 


As previously stated, the development of kinetic mechanisms 
should be considered from the differential reaction rate viewpoint 
rather than from the steady-state equilibrium ratio viewpoint so 
prevalent in the standard chemical engineering texts on the 
subject. Only in this way may the true relationships of some of 
the transitory species such as radicals be properly considered. 
The previous reference“ has described more completely the 
procedure which should be followed in determining the various 
steps in a kinetic mechanism by analog simulation. 


Even with the help of the analog computer this is still a very 
tedious task since a complex trial and error procedure must 
be used to determine the best values of the reaction rate con- 
stants. Therefore, the recent work by the International Business 
Machines Corporation on digital computers and by Computer 
Systems Incorporated“? on analog computers promises to 
eventually make possible the determination of the proper 
values of kinetic constants automatically, by a machine con- 
trolled and operated trial and error procedure. The investigator 
can then try a very large number of possible mechanisms in a 
short period of time. Since the machine also determines the 
value of the minimum root mean square error function for each 
separate mechanism, it is a simple matter to then choose the 
most probable mechanism from the available choices presented 
to the computer. 


Batch reactors and stirred tank continuous reactors along 
with their kinetics are readily simulated on large scale analog 
computers including | their control relationships provided only 
that reactor mixing is reasonably complete. However tubular 
reactors are another matter if their control functions must be 
considered. Because of the multiple sectioning necessary for a 
finite difference type approximation in real time, most computers 
are quickly overtaxed. 

The only solution is then to use a repeated multi-step 
procedure. The steps required are: 


(1) Determine the steady-state composition distribution in 
the reactor by converting the length variable to time and solving 
on analog or digital computers. 

(2) From the resulting composition history determine the 
regions of possible heat exchange and the heating or cooling 
load required. 

(3) Solve the reactor dynamically (i.e., time-wise) for the 
heat transfer-control relationships only. If the control system 
cannot regulate the required heat transfer from step 2 as an 
input function at the Proper point determine the resulting 
temperature excursions in the reactor. 

(4) Repeat step 1 for the limits of the temperature excursions 
found in step 3. Determine whether or not this causes the outlet 
composition to deviate outside the allowed limits. 

(5) If the composition deviation is too great determine the 
new heat load requirements for the temperature control system. 

(6) Revise the control system as necessary and then repeat 
step 3 to determine new temperature excursions. 
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(7) Repeat as necessary until the temperature contro] system 
can properly regulate the outlet composition for the reactor. 

While this is a difficult, long-term procedure it is the only 
one known which can be used at present to check the suitability 
of a tubular reactor control scheme by theoretical means. 


Other Chemical Engineering Computations 


Just as was the case with reaction kinetics and reactor 
dynamic investigations, most other types of chemical engineering 
calculations are still in a state of flux. For example, despite 
years of effort and unnumerable papers on the subject, the 
design of distillation columns by trial and error algebraic 
solution of the steady-state equations of the column is still 
unsuccessful in a large proportion of the important cases. This 
is particularly true for those cases where multiple feeds or 
sidestreams are present, where very close separations are 
necessary, or small impurities are present, or where non-ideal 
equilibrium relations are present. Unfortunately these are the 
cases of greatest interest and at the same time the most difficult 
to carry out by hand methods. 

There is however a method of machine computation, which, 
though tedious, can be used for such computations. This is an 
integration of the dynamic equations of the column, a set of 
ordinary differential equations“. The steady-state condition is, 
of course, a special case of the solution of the dynamic behavior 
of the column. 

Computation of the design of absorption and extraction 
columns is an even more difficult task than that of distillation 
columns since the equilibrium relations involved are still more 
non-ideal than in distillation. Again, though, the integration 
scheme will work if the equilibrium relations involved are known. 


THE TOOLS AND THE PEOPLE REQUIRED 


With some knowledge of the places where computers can 
and are being applied one may discuss the required computer 
sizes for specific tasks and the numbers and types of personnel 
required to operate them effectively. 


Some Notes Concerning the Required Size of an 
Analog Computer Installation 

There is a fact of life about all computers, particularly 
analog computers, which should be kept in mind by everyone 
concerned with their procurement and use. This is the fact 
that no computer (at least to date) has even been big enough to 
handle all of the problems which would eventually be applied to 
it in the future. As a result, computer facilities have continued 
to grow and grow in size, flexibility, and speed in a.most a 
geometrical progression. 

Therefore, the first caution one can give to the organization 
contemplating a computer installation is to be sure that the 
resulting machine and installation are readily expandable when 
the need arises. 

There are four main uses to which a chemical industry 
analog computer could be placed. These are: 

(a) Instruction in computer usage, process dynamics, and 
control system operation. 

(b) Investigations of control system (i.e., instrumentation 
only) response and design. 

(c) Simulation of control system response and design 


‘(including rudimentary process dynamics and chemical kinetics). 


(d) Complete simulations of chemical kinetic systems, 
process development investigations, and overall processing 
system designs with control included. 

Each of these in turn requires a larger and much more 
versatile computer as shown in Table 3. Of course, the larger 
machines can readily carry out the tasks of the smaller machines 
and separate installations are not required, unless the workload 
of the larger machine should necessitate such action. 

The first thing the reader will note concerning the computers 
recommended in Table 3, particularly if his experience has been 
with computers in the defense industries, will be the extremely 
high proportion of non-linear gear, especially in the larger 
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TABLE 3 


DISTRIBUTION OF COMPONENTS FOR TYPICAL 
ANALOG COMPUTER APPLICATIONS IN THE CHEMICAL INDUSTRY 
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Operational Amplifiers | 10-20 | 15-50 | 40-100 | 80-200+ 
Multiplication Channels | a 2-5 10-30 | 50-100+ 
Function Generation Channels | 1-2 2-6 5-10 | 25-50+ 
Operational Relays |} 1-2 6 12 | 25-50+ 
Function Switches | 3 6 12 | 25-SO+ 
Separate Control Consoles | 1 1 1-2 |} 2-3+ 
Plotting Channels x, y | 1 | 1 | 3 3-4+ 
Plotting Channels x, t 2-4 6 12 | 18 
Repetitive Operation Yes | If Desired | If Desired | If Desired 
| | 
machines. These recommendations are made because of the 


nature of chemical problems, particularly those involving chemi- 
cal kinetics. Here a very large number of multiplications and 
function generations, particularly exponential functions, are 
required. In addition, the large selection of operational relays 
and function switches permits the direct simulation of multi-step 
processes — such as batch drying and crystallizations as a 
continuous operation — thus avoiding the necessity of breaking 
the problem into a number of separate steps w ith the attendant 
possibilities of error and difficulties of correlation of results. 


Other points which should be particularly emphasized con- 
cerning Table 3 are as follows: First, repetitive operation is 
extremely useful for instructional purposes, since its use with 
an oscilloscope display permits the most extensive display of 
the dynamics of a process or control system in the shortest 
possible time. In addition, this “continuous display”’ demon- 
strates most vividly to the student the effect of changes in the 
sy stem parameters and thus gives him by far the best intuitive 

“feel” for process behavior. However, there are severe limita- 
tions on plotting accuracy with such displays. Thus, they may 
be used for research and design problems only as rough ‘ ‘ball- 
park” approximations. The solution should be repeated at 
normal speeds to get the final most accurate answer. Thus, 
repetitive operation is not specifically recommended for the 
larger installations but is left a matter of personal preference 
of those involved. 


Second, the extreme flexibility of the x, y plotter as a function 
generator should be brought to the readers’ attention. It is 
particularly useful for presenting such functions as temperature 
variable heat transfer coefficients which cannot readily be pre- 
sented by standard function generation equipment. 


Organization and Numbers of Personnel Required 


Just as each of the different contemplated applications of 
the analog computer requires a somewhat different size of 
computer and a different distribution of components, a some- 
what different organization is necessary for their proper utiliza- 
tion. Table 4 outlines recommendations in this regard. 


Several points concerning Table 4 are worth reiteration or 
further explanation. First, since electrical and mechanical 
engineers at present have a better background in automatic 
control theory and the associated mathematics, they are 
recommended almost exclusively for the smaller organizations. 
However, as the group delves more and more into the i intricacies 
of the chemical processes themselves, the chemical engineer's 
knowledge becomes more and more important until they out- 
number the others. Note that a certain number of electrical 
engineers are still recommended because of their familiarity 
with the operation and capability of the computer. 


Second, full-time maintenance personnel are recommended 
only when the computer becomes large enough to demand such 
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TABLE 4 


NUMBERS AND BACKGROUND OF PERSONNEL REQUIRED 
To OPERATE _ANAI ALOG / Compu TER FACILITIES OF VARIOUS SIZES AND USAGES 
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Analog Computer Technicians 7 — 1-2 | 2-3+ 
(Repairmen) 
Open Shop |} Yes | Yes | As Desired | No 
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attention (50-100 amplifiers). For smaller machines, the attend- 
ant engineers should be capable of performing the necessary 
reventive maintenance checks and repairs. 
Third, the staff for the larger computer continues to increase 
even though the contemplated problem size has also increased. 


This is due to the fact that the economics of the utilization of 


large-scale computers almost demands their use on a multi-shift 
basis because of the very high depreciation charges involved. 
Thus, the staff recommended here is for two shift operation of 
the largest machine recommended. 

Fourth, because of the difficulty of training the necessary 
personnel and the extreme complexity of the problems involved, 
it is the majority opinion of those concerned that a very large 
computer facility should be operated on a closed- -shop basis. 
However, every effort should be made to obtain the complete 
cooperation of the person originating each of the problems 
solved, even to having him observe the actual solution on the 
computer. Often he will have a “feel” for the problem which 
will be invaluable in checking out the computer and evaluating 
the final results. 


New Developments in Analog Computers 


There are several recent or eminent developments in analog 
computer equipment and in operating techniques which promise 
to greatly aid and possibly even revolutionize analog computer 
operation. Most of these involve digital techniques of one type 
or another. 

One of the major present-day problems with analog computer 
operation in the process industries is the difficulty of adequately 
representing a true time delay. The dev elopment of a suitable 
and reasonably priced magnetic tape or drum mechanism for 
this purpose would be a major boon to the field. Such devices are 
now under development. 

Analog computers are especially designed for the solution 
of ordinary differential equations and do this extremely well. 
Digital computers on the other hand are very fast and accurate 
in the solution of algebraic equations but relativ ely slow with 
large systems of differential equations. In addition, analogs 
have only rudimentary logic capabilities, while extremely com- 
plex logic operations are carried out directly on digital machines. 
There is thus a very fruitful potential in the joint use of these 
machines for the solution of large process development-type 
problems which entail all three types of mathematical operation, 
solution of differential equations, solution of algebraic equations, 
and logical selection from several possible paths of action. 
Devices for connecting the two computers together during 
operation are in use by several of the defense system contractors, 
such as Space Technology Laboratories, Convair Astronautics, 
and International Business Machines. Models of these are also 
under development for use by the process industries. 

Another area where aid is needed in analog computation is 
in the simulation of distributed parameter systems and processes. 
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Adequate mathematical representation of these systems requires 
the use of partial differential equations or extremely large sets 
of simultaneous ordinary differential equations. Thus, their 
solution on analog equipment requires a very large machine. 
However, the use of a modified form of the data storage drum 
mentioned above for true time delay simulation offers an aid in 
the solution of this type of problem. If preliminary results can 
be stored during the solution of a problem of this type on the 
analog computer, it can be solved on a much smaller machine 
than required above by the use of a repetitive trial and error 
procedure which can be readily programmed. Again these 
devices are under active dev elopment by the computer manu- 
facturers. Thus, digital devices promise to have a major effect 
on analog computer operation in the process industries. 


Digital Computer Installations 


One should not leave the discussion of computers, particularly 
after presenting details of analog computer sizes and the per- 
sonnel required, without presenting similar data concerning 
digital computer installations. 

Digital computers, because they are available in only a few 
specific models, do not present as difficult a choice in determining 
the best size of installation for a specific purpose as do analog 
computers, since one need not break them down to specific 
components. The best rule here is to procure the biggest and 
fastest machine that one can finance and staff. Just as in the 
case of analogs, the company invariably underestimates the 
applications that will be made of the computer once it has been 
installed. 

Computer installations of the type necessary for calculations, 
such as those described in this paper, can be divided into two 
categories, depending on the size computer installed. 

These are first, the medium sized computers of which the 
Royal McBee LGP-30 is considered the minimum size. These 
progress through the Bendix G-15, and the IBM-650, to the 
Burroughs-Electrodata 205. The second group, large sized 
computers, starts with the Burroughs-Electrodata 220 and pro- 
gresses through the IBM 704, 709, and Remington Rand Univac 
1103A, to the RCA 501, and the recently announced IBM 7090. 

Table 5 presents the author’s opinions as to the personnel 
required to adequately man a digital computer of either size. 
The ranges in numbers of personnel recommended are to allow 
for the variations in computer size within each group. 


TABLE 5 


PERSONNEL REQUIREMENTS FOR TYPICAL PROCESSING 
Ibu STRY DiGITaL COMPUTER INSTAL LATIONS 











Large Size 
Burroughs 
Datatron 220 


Medium Size 1BM 
Royal McBee 704 
Size of Computer LGP-30 Remington Rand 
Involved Bendix Univac 1103A 
—_—_— G-15 M 
Type of Personne! IBM 709— 
Considered 650 Remington Rand 


Burroughs 


Univac 1105 


Datatron-205 IBM 
7090 
RCA 
501 
Statistical Specialists 
1. B.S. Ch.E. with 
M.S. in Statistics 
1-2 2-3 


or 
2. Ph.D. Ch.E. with 
Statistics Minor 
Chemical Engineering 
Design Specialists 
M.S. or Ph.D. Ch.E. | 2-3 | 4-8 
with Applied | 
Mathematics Minor 
Linear Programming: and | 
Operations Research 
Specialists | 
1. M.S. or Ph.D. Ch.E. | | 
| 
| 
| 





with Applied 0-2 3-8 
Mathematics Minor | 
2. M.S. or Ph.D. | 
Mathematiciaas | 
with B.S. in Ch.E. | 
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THE IMPACT OF COMPUTER CONTROL AND A DIS- 
CUSSION OF ITS PROPER PLACE IN THE CHEMICAL 
INDUSTRY 


Present Status 


Digital computer process control is, as everyone knows, 
one of the most discussed, evaluated, and predicted subjects 1 in 
the industry today. Therefore it deserves a special mention and 
evaluation here. Nearly every meeting of any one of the various 
societies has its full share of such papers — most of them of a 
speculative nature. Despite this, surprisingly little has yet been 


published or presented concerning the true performances of 


these machines in chemical plants. In addition, because of digital 
computer control, the past two and a half years has been a period 
of some of the greatest duplication of effort on the part of the 
computer manufacturers and of the processing companies in our 
history. Each of the approximately 16 computer manufacturers 
who build these machines and the 20 processing companies who 
have purchased them has — by means of both internal and inter- 
industry study projects — sought out for itself where the 
economic returns from these new devices may lie. Proprietary 
consideration, necessary on both sides to protect any and all 
economic advantages uncovered, has prevented and still does 
prevent a complete and free disclosure of the results of these 
tests. Only in very recent months have articles devoted to 
specific digital computer installations begun to appear“ *:” and 
even these must be somewhat less than full discussions of the 
topic, for the same reason. 

There is, of course, no remedy possible or even desired to 
such a situation since our economic system demands it. It is 
mentioned here only because of the drain it has made upon the 
very small pool of trained manpower available for all process 
control and process dynamics work. 

In this discussion no attempt will be made to compare the 
specifications and characteristics of the various models since this 
has been admirably handled by Freilich for all but the very 
newest of the machines. 

Much of the available information concerning process control 
computers and their applications can be summarized as in Table 6. 
Since these machines have been used exclusively for supervisory 
control functions, their place in chemical process control can be 
diagrammed as shown in Figure 1. Along with the job of main- 
taining an economic balance of the plant’ S operation, there are 
several additional tasks to which they may readily be applied 
without hindering the performance of their primary function of 
optimization. These are listed in Table 7. 

It is possible now to develop some basis for economic and 
technical justifications which can be put forward for the use of 
process control computers as steady state economic optimizers 
in chemical processing plants and petroleum refineries. 

Where are Computers Justifiable— 
Economic Considerations 

If one can make an assumption of generality concerning 
some conclusions that have been drawn concerning the applica- 
bility of process control computers, one can readily diagram 
the relationship of digital process control computer installed 
cost to the process plant total investment and the required 
return on total investment before taxes which is attributable 
directly to the work of the digital control computer. 

The resulting general statements are: 

(1) A return on total plant investment of at least 0.5% per 
annum is probably possible on any process due entirely to the 
monitoring and long term economic optimizing a abilities of the 
computer alone. 

(2) If a plant systems study shows that a return of greater 
than about 6% per annum is possible, serious engineering study 
should be given to the possibilities of removing process bottle- 
necks or improving process work in the plant. After this, the 
plant should be restudied to determine the remaining return on 
investment which can still be credited to the computer. If such 
returns are possible the plant was obviously not designed origin- 
ally for optimum compatibility of its constituent units. 
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TABLE 6 
SUMMARY OF AVAILABLE INFORMATION CONCERNING 
Process CONTROL DiGITAL COMPUTERS 


(1) All installations to date have exercised supervisory control of the process 
only, either on-line or off-line. That is, they have been used to adjust the 
set points of conventional instruments either directly from the computer 
or through manual intervention. 


(2) All computers so far have been used for steady-state analyses of the process’ 
performance, mainly analyses of the economic returns possible from the 
process. Thus these computers have not been concerned with the stability 
of the >rocess control itself or with process start-ups or shut-downs 


(3) Being mainly solid state and derived from tested weapons system or large 
commercial computers, the process control computers are showing as a 
rule a truly remarkable reliability and dependability. 


Study is showing that the earlier machines contained too small a memory 

capacity to be truly satisfactory process control computers. Thus the 
more recent models have provision for large and quite flexible memories. 
Also several of the older models are being redesigned or modified to have 
this feature added. 


(4 


~ 


TABLE 7 


TASKS FOR WuicH DiGitaL COMPUTER CONTROL MAY BE REQUIRED 
Prov IDED  E CONOMIC Ju STIFICATION IS 3 POSSIBLE 


(1) Modification of plant mathematical models on a regular basis to account 
for such effects as rapid catalyst deterioration, tar buildup, heat exchanger, 
condenser, or column fouling and others of a similar nature. These coeffi- 
cients in the plant economic control model must change significantly over 
periods of as short as a few or even one day to preclude the use of off-line 


computers for such tasks. 

Where the plant mathematical model is too complex, i.e., the plant inter- 
actions are too complex, or their kinetics are too detailed, etc., for special 
analog type computers to be designed for the task. 

Where the process or kinetic model is incomplete in that statistical analysis 
must be regularly made to determine and upgrade kinetic and other 
mathematical model coefficients and terms. 

Data logging, data monitoring and data processing and evaluation. Such 
tasks are subsidiary to the digital computers main task but may add 
justification for its use. 

Statistical analysis of plant operating data where these must be carried 
out on a regular repetitive basis. 

Regulation of cyclic operations where variable time periods within the 
cycle are necessary. Batch operations employing catalysts of varying 
activities are of this type. 


Emergency procedures such as partial plant shutdowns. 
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While the actual limits used here may vary from company 
to company depending upon accounting and management prac- 
tices, a relationship such as suggested above appears to hold 
throughout our industry. W ith it in mind, Figure 2 outlines 
the allowable computer size versus plant investment if only the 
computers monitoring ability is to be capitalized on and none 
of its special control abilities as listed in Table 7. Lines are 
shown for one, two and three year payouts after taxes. 

Figure 3 details the range of allowable computer installed 
costs versus plant total investment between the limits of 0.5 
and 6% return, with the locations of intermediate return rates 
indicated. This figure is for a three year payout. Finally 
Figure 4 lays out in regions of percent return versus plant 
overall investment the areas where computers are readily justi- 
fied, where they must be justified by specific studies, where 
further plant studies for removing bottlenecks should be made, 
and where computers can not be justified at all because plant 
investment is too low. Figure 4 shows what is probably the 
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Figure 1—Present day concept of the place of computer 
control. 
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Figure 2—Allowable total cost of computer installation to 

capitalize on the computers monitoring ability alone 

(assuming 0.5% return on original investment due to 
monitoring ability of computer ). 


present day standard installation, that is, one costing about 
$300,000 complete and requiring a two year payout after taxes. 
I q $ y ; 


The most obvious points brought out by these graphs are 
first, the extreme dependence upon plant size or better plant 
total investment, and second, the inverse effect of computer 
installed cost. From this one can immediately conclude two 
things; that computer control will always be more important to 
the petroleum refineries than to chemical plants and that computer 
control will be much easier to justify on new plants. The first 
is true because of the generally much larger size of the refineries. 
Also within the chemical companies, petro-chemical installations 
will be the most susceptible, again because they are usually 
larger. New plants are more susceptible to computer control 
than old because here one is not faced with a tremendous outlay 
of capital for instrumentation just to make the computer and 
plant compatible. If properly planned ahead of time, the 
instrumentation for a new plant to allow computer control 
should be little more expensive than that planned to operate 
without it. Thus in contrast to the older plant (where the capital 
investment includes both computer and instrumentation), the 
required investment in the new plant includes only the computer 
and a small amount of associated hardware, and may be only 
one-third to one-half that of the former. 
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There are also several other factors which should be kept 
in mind while establishing the economic justification for computer 
control. These are listed in Table 8. 


Where are Computers Justified— 
Technical Considerations 


Considered along with the economic justification of digital 
co:nputer control and directly related to it is the equally impor- 
tant question of technical justification, i.e., is the computer 
really necessary for the task assigned to it or will a simpler, 
cheaper control method accomplish the same results. ‘Table 7 
has already listed several tasks which the digital computer can 
perform and for which it may be necessary. Some other 
considerations which should be kept in mind are listed in Table 9. 

This is also the place to emphasize the necessity for making 
a thorough plant systems study and for developing a mathe- 
matical model of the plant which is as complete and as exact 
as possible”, This is necessary both to permit the proper 
programming of the machine and also to thoroughly document 
the claimed economic justification for its use. 

It has often been remarked by systems people that there 
has never been a thorough plant study which has not paid for 
itself in new knowledge concerning the process and its operation 
whether a computer was actually purchased or not. While only 
about one-fourth to one-fifth of the chemical plant studies to 
date have uncovered potentially profitable digital computer 
applications“ their uncovering of new plant knowledge has 
made them a very worthwhile endeavor for the systems groups. 


TABLE 8 


ADDITIONAL ECONOMIC FACTORS IN 
JustiFYING DiGiITAL COMPUTER CONTROL 


(1) One of the primary findings of most systems studies of existing processes 
is the great gains which can be made from the installation of analytical 
instrumentation alone. Since such instrumentation is usually necessary 
to allow digital computer control to even be applied, economic gains really 
attributable specifically to the more modern instrumentation are often 
credited to the computer. Great effort should therefore be made to separate 
those gains due mainly to the greater process knowledge available from 
analysis instrumentation from that which is due strictly to the computers 
use of this information. 


Similarly one should be careful that the computer installation is not credited 
with gains attributable only to the specific new knowledge and resulting 
adjustment of plant operating methods which frequently occurs as a result 
of a complete plant systems study and which would be applied whether 
the computer had been purchased or not. 

(3) For cases where some of the coefficients of the plant mathematical model 
are unknown or variable and it is contemplated that the computer will 
be necessary to compute them, consideration should be given to the situation 
which will exist once the computer has determined these data or established 
a correlation for it. Could the plant be converted to some simpler control 
system at that time? If so, serious consideration should be given to 4 
temporary installation of the computer using one of the several rent 
systems available rather than an outright purchase of the machine. 
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TABLE 9 


SoME TECHNICAL CRITERIA FOR 
JUSTIFYING COMPUTER Con TROL 


The computations which the computer is to carry out must be of such complexity 
that no simpler device could perform them sufficiently well to capture all or 
most of the contemplated returns. That is, the following possibilities must be 
considered and discarded as impractical for the reasons named: 
(1) Nomographs — Solution of simple algebraic equations can often be carried 
out by means of nomographs or tables. Thus a certain definite level of 
complexity is necessary to justify the computer. 
Analog Computers — Systems of differential equations or of algebraic 
equations expressed in differential form are rapidly solved on analog 
computers. For the type of control tasks contemplated today this way 
will be much cheaper than the corresponding digital machines. However 
if the optimizing routine includes the adjustment of several of the coefficients 
of the plant mathematical model, analogs will be impractical and the digital 
machine will be necessary. 

(3) Off-Line Digital Computers — If the computatio.as carried out by the 
control computer must be repeated at intervals of less than about one 
hour, off-line computation would be impractical, otherwise it must be 
considered a definite possibility and a factor in the justification of the 
on-line, in-plant machine. 
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It has been mentioned by some in taking a quick survey of 


the field of computer control that the chances of justifying a 
computer are inversely proportional to one’s degree of knowledge 
of the process and its kinetics. While this is of course an ex- 
aggeration it is nevertheless true that the necessity for making 
a statistical analysis as part of the computer’s task in optimizing 
the process greatly increases the possibility of a justification 
for its use. 

If off-line analog or digital computation can uncover — by 
means of theoretical considerations — the kinetic or mass trans- 
fer relationship which is the basis for the required statistical 
evaluation, then a simpler control system is often entirely 
sufficient. 


The Place of Analog Computers in Process Control 

Present day process controllers can be considered as essen- 
tially one input and one output systems. Therefore one can 
only control a multi-variable system with conventional equipment 
through the use of multiple controllers in a cascade arrangement. 
However, special analog computers can now be designed to 
carry out these and even more complex control functions in one 
“black box”. Two of these — one designed to correct for the 
effect of ambient conditions on the operation of a fin-fan cooler 
used as a condenser on a distillation column”, and the other 
designed to make an overall heat balance about a polymerization 
reactor“) — have already been described in the literature. 
The former is being manufactured by several instrument and 
computer companies. The author’s investigations of processes 
have to date revealed many potential applications for such 
specialized controllers, particularly in the control of chemical 
reactors and in the modification, interpretation, and presentation 
of plant data. Also in contrast to the larger digital computer 
installations (which because of their complexity and size are 
usually one of a kind projects), the analog computer-controller 
designs are often directly transferrable from one process to the 
next. Because of the press of their other duties, most of the 
systems and control groups have not had the time, manpower 
or the facilities to develop these items themselves. It is therefore 
hoped that the instrument and computer companies will eventu- 
ally take a more active interest in this presently somewhat 
neglected field. It is this author’s feeling that this can be a 
very active sales field in the near future once study has shown 
the best form which these devices should take for maximum 
flexibility of application from one process to another. 


Dynamic or Real-Time Control by Digital Computer 


So far digital computer control has only been applied to 
static or steady-state control of chemical processes. What about 
its use for controlling processes which are subject to rapidly 
changing operating conditions such as cyclic or batch reactors 
or tubular reactors with variable feed conditions? The reader 
should particularly note that the examples considered in this 
question have all been reactor situations. As will be brought 
out later, only reactors appear to have the non-monotonic or 
exponentially monotonic behavior which could possibly necessi- 
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tate in-line dynamic control. All other processing operations 
appear to exhibit an overdamped monotonic response which 
can be adequately handled by conventional instrumentation with 
the possible aid of supervisory optimizers. Note that reactor, 
as used here, is the general term which includes the reforming 
and cracking, etc., operations of the petroleum refinery as well 
as the various reactors and furnaces in the chemical piant. 

In an attempt to evaluate the above question, a generalized 
chemical processing problem has been formulated and presented 
to the computer manufacturers for their use“. Study of this 
model has revealed that dynamic control of such situations with 
the inclusion of an economic optimization feature will be an 
extremely complex calculation even for relatively simple 
cases5.16.12, Thus present day process control computers are 
probably much too small (i.e., too little storage) and too slow 
in their computing speed to be used for this purpose. 

At the same time it must be kept in mind that this is the same 
area for which several possible analog computer controllers 
have been postulated. ‘Thus there is the possibility that an 
analog computer controller added to the reactor’s internal 
control loop may make the reaction process perfectly amenable 
to overall supervisory digital computer control of an economic 
optimizing type, if the digital computer is used to reset some of 
the potentially variable control functions of the analog computer 
controller. Because of its dynamic control ability, ‘the analog 
computer controller may readily handle the reaction process 
while the digital computer calculates the next group of optimum 
settings at its leisure. 

It is, of course, also possible that continued research into 
dynamic optimization procedures and in kinetics and reactor 
dynamics will reveal an economic justification for a larger and 
more capable digital computer for these specific process ‘control 
duties. In this case one of the several new solid-state scientific 
machines may serve admirably and preclude the necessity for 
the development of a larger special control computer. 


WHAT THE CHEMICAL ENGINEER MUST DO 


Because of the sudden breaking of the field of computer 
analysis and its rapid dev elopment of the frontiers of science 
combined with the chemical engineer’s previous devotion to 
the steady-state aspects of process plant design and plant 
technology, he has found himself ill prepared to capitalize on 
computers and their capabilities. This combined with a differ- 
ence in industrial emphasis has resulted in a decided difference 
in the degree of application of computers in our branch of 
engineering compared to aeronautical, electrical, and mechanical 
engineering. 

To make the best use of computers, the chemical engineer 
must have a broad knowledge of basic chemistry, physics, and 
properties of materials. He must be thoroughly grounded in 
mathematics with a good facility in one or preferably more of 
the following: operational mathematics, numerical analysis, and 
statistics. In addition, he must be able to think in terms of 
dynamic behavior of chemical processes: of reaction kinetics, 
of fluid dynamics, and of mass and heat transfer. 

These needs are resulting in a rapid revamping of the chemical 
engineering curricula to provide a much more fundamental 
program w ith more emphasis on mathematics, chemistry, 
physics, and control system response. To aid the engineer 
who is already in industry, the technical journals have published 
many tutorial series of articles and bibliographies“ on the new 
aspects of chemical engineering. In addition, the engineering 
societies are taking an active hand in supplying such information. 
Examples of these activities are the Summer Conferences on 
Control Theory for Chemical Engineering Teachers and the 
Special Lecture Series of the American Institute of Chemical 
Engineers. 

In addition, several of the computer manufacturers such as 
Electronic Associates Inc., and the International Business 
Machines Corporation (IBM) conduct periodic short courses 
in computer programming and computer application to process 


225 





industry problems. Also it is the practice of practically all 
computer manufacturers to train a certain number of the cus- 
tomer’s personnel to operate and maintain each new computer 
purchased. Thus there is every opportunity for the chemical 
engineer who so desires, w hether still a student or already in 
the field, to acquire the necessary knowledge and experience 
to take his place in the rapidly growing field of computer usage 
in our industry. 


WHAT ABOUT THE FUTURE? 


Approximately three years ago this author made a prediction 
that the year 1968 would see the followi ing in the use of computers 
in the process industries: 

(1) 75 of the 100 largest process companies in the United 
States would be engaged in simulation studies with analog 
computers or advanced digital computers. 

(2) All Jarge chemical and petroleum companies would have 
large scale digital equipment devoted to scientific and engineering 
work. Many of them will have several of these machines. 

(3) If instrumentation development keeps pace, 150 or more 
large scale machines will be devoted to repeated plant optimiza- 
tion studies. Many will be connected as overall plant supervision 
machines. 

(4) Industry will be actively engaged in ‘‘economic automa- 
tion” or the maximum of automatic control which will show an 
economic payout. 

With the ten-year period only a quarter over these predictions 
appear, if any thing, conservative and we are well on the way to 
fulfilling all of them. This is particularly true of the last two 
since instrumentation developments have also been quite spec- 
tacular in the last two or three years and promise to continue. 


References 

(1) Williams, T. J., “Systems Engineering for the Process 
Industries—Its Application, and Its Future Potential”, 
Fourth E. P. Schoch Lectures, U. of Texas, Austin, Tex., 
October 16-17, 1959, currently being published in serial 
form by Chem. Eng. 


226 


(2) 


(3) 
(4) 
(5) 
(6) 
(7) 


(8) 


~— 


(9) 


(10) 
(11) 


(13) 


(14) 


(15) 


(16) 


(17) 


(18) 


Machine Computations Committee, Am. Inst. Chem. Eng! 
“Computer Procedures for Estimating Physical Properties”, 
September 1960. 

Petersen, T. I., Chem. Eng. Progr. Symposium Ser., 56, 
No. 31, 111-120 (1960). 

Computer Systems Inc., Form 80-105-003, 
Junction, N.J., 1960. 

Greenstadt, J., Bard, J., and Morse, B., Ind. Eng. Chem, 
50, 1644-1647 (1958). 

Rosenbrock, H. H., Brit. Chem. Eng., 3, 364-367, 432-435, 
491-494 (1958). 

Eisenhardt, R. D., and Williams, T. J., Contr. Eng., 7, 
No. 11, 103-114 (November 1960). 

Madigan, J. M., Chem. Eng. Progr., 56, No. 5, 63-67 
(May 1960). 

Parsons, J. R., Tolin, E. D., and Andrews, A. J., Contr, 
Eng., 7, No. 9, 150-153 (September 1960). 

Freilich, Arthur, ISA Journal, 6, No. 7, 47-65 (July 1959). 
Hall, C. R., Chem. Eng. Progr., 56, No. 2, 62-66, (Febru- 
ary 1960). 

Lupfer, D. E., and Berger, D. E., ISA Journal, 6, No. 6, 
34-39 (June 1959). 

Tolin, E. D., and Fluegel, D. A., ISA Journal, 6, No. 10, 
32-38 (October 1959). 

Williams, T. J., and Otto, R. E., Am. Inst. Elec. Eng. 
Paper CP-60-119, Winter General Meeting, Am. Inst. Elec, 
Eng. February, 1960. 

Beecher, A. E., and Gould, L. A., Am. Inst. Elec. Eng, 
Paper CP-60-980, Joint Automatic Control Conference, 
Cambridge, Mass., September 1960. 

Boydston, R. E., Am. Inst. Elec. Eng., Paper CP-60-979, 
Joint Automatic Control Conference, Cambridge, Mass., 
September 1960. 

Kipiniak, W., and Gould, L. A., Am. Inst. Elec. Eng,, 
Paper CP-60-981, Joint Automatic Control Conference, 
Cambridge, Mass., September 1960. 

Williams, T. J., Chem. Eng., 68, No. 1, 103-110, January 
9, 1961. 


Monmouth 


x &k * 


The Canadian Journal of Chemical Engineering, October, 1961 








A 
Gener 
leasec 
for tl 
shipps 
89,001 

Ne 
ration 
ata sz 
to the 


Comn 


he 
ar 
pressur 
tremeni 
effect, 
As 
commo 
conditic 
highly | 
water | 


Econo 


In t 
acid sol 
pipe lit 
econom 
to prod 

Bec 
the size 
of shipy 
in sing] 
weights 


Design 

E'xa 
clearly 
and at : 
By com 
of 166 | 
low ter 
commox 
econom 
ence on 


1Manusce1 
2Presiden 
Montrea 


The Ca 





” 


-” 


g. 
Cc. 





61 


Rail Transportation of Hydrogen 
Chloride Anhydrous’ 


Cc. H. DRURY? 


A railway tank car has been built for Canadian 
General Transit Co. Limited and is operating in the 
leased service of Dow Chemical of Canada Limited 
for the shipment of hydrogen chloride anhydrous 


shipped at minus 50°F. in single shipments of about 
89,000 Ib. 


Necessary safety devices and a liquid CO, refrige- 
ration system are included which will maintain the car 
at a safe pressure level for an extended period subject 
to the supply of liquid CO.,. 


Commodity 


Zs principal detailed phy sical data of hydrogen chloride 
anhydrous are shown in Figure | together with its vapor 
pressure-temperature curve. In addition, this commodity has a 
tremendous affinity for water and, because of its dehydratin 1g 
effect, is dangerous to living tissue. : 

\s a dry gas it is essentially inert and does not attack the 
commonly-used structural metals under temperature and pressure 
conditions normally used. In the presence of moisture it becomes 
highly reactive w hich calls for special design where contact with 


water vapor is possible, v valve stems, open connections, etc. 


Economic Considerations 


In the past, the cost “of transporting aqueous hydrochloric 
acid solution combined with the fact that no means, other than 
pipe line, was available to transport the anhydrous materiaj 
economically forced many users of HCI to locate plants adjacent 
to producers of HCl. 

Because unit costs of shipment normally vary inversely to 
the size of the shipment it was decided to inv estigate methods 
of shipping substantial weights of hydrogen chloride anhydrous 
in single shipments by rail and in order to obtain economic 
weights i it was obviously essential to ship in liquid form. 


Design Considerations 

Examination of the HC] temperature vapor pressure curve 
clearly indicates at 90°F. a vapor pressure of 795 p.s.i. gauge 
and at about minus 45°F. a vapor pressure of 100 p.s.i. gauge. 
By comparison, anhydrous ammonia at 90°F. has a vapor pressure 
of 166 p.s.i. gauge. With such characteristics it is evident that a 
low temperature technique of handling was essential. If the 
commodity can be held at that temperature the problem of an 
economic railway shipping container is relatively simple. Experi- 
ence on low temperature rail shipments of liquid carbon dioxide 


1Manuscript received April 12; accepted July 12, 1961. 
2President and Managing Divwoter, Canadian General Transit Co. Limited, 
Montreal, Que. 
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was available and it appeared that some adaptation of this 
container would be satisfactory. 

The railway tank cars approved for handling HCI anhydrous 
are constructed of steel to Specification A-300, Class 1, ASTM 
\212, Grade B, suitable for service to minus 50°F. and are 
specially and carefully insulated with 10” Sty rofoam, the whole 
enclosed in a steel jacket, the entire design requirements coming 
under Container Specification ICC-105 A-600W. (Hy draulically 
tested to 600 p.s.i. with safety valve settings at 450 p.s.1.). 

In order to ensure no liquid will be discharged through the 
safety valve the ICC limits the allowable filling density (ratio 
of the weight of commodity in pounds to the water capacity of 
the container in pounds) for HCl to 85%. (ICC Special Permit 
No. 2308; BTC Special Permit No. 467). Using this factor 
and a net cubic capacity of about 10,500 Imperial gallons, 
approximately 89,000 pounds of hydrogen chloride can be carried. 


Heat Balances 


In order to establish how long the car would remain safely 
in an ambient temperature of 90°F., on the assumption that the 
lading was shipped at minus 50°F., a simple heat balance was 
done which indicated a time lapse of minimum 30 days before 
the vapor pressure of the hydrogen chloride would reach the 
discharge pressure setting of the safety valve at 450 p-S.., 
equivalent to a temperature of approximately 47°F, 


1200 HCL T2124 7°F 
Pe #185 psig 


CO,> & = 87. OF. 
R #1051.6 psig 
Hoo L 


1000r 
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700+ MELTING POINT -101%C (-168°F) 

DENSITY @0°C.760 mm. Hg-gm/iter!.639 o 
600. HEAT OF VAPORIZATION@ -85°C. cal/gm™m 105.9 v 


soo- & 
400- & 


300+ 


PRESSURE-IN POUNDS PER SQUARE INCH GAUGE 


200- 
100} - 
L 1 1 1 —a 1 1 L L L L 1 1] 
-140 -120 -100 -60 -60 -40 -20 t0 20 40 60 eo 100 120 
TEMPERATURE -IN DEGREES FAHRENHEIT 


Figure 1—Properties of hydrogen chloride, anhydrous. 
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Safety 

Although 30 days would seem reasonably adequate protection 
it was decided that under extreme conditions this time might be 
exceeded and therefore it was desirable to provide additional 
protection by refrigeration which could maintain the lading at 
a suitable temperature for an additional period. 

A number of systems, including mechanical refrigeration, 
were considered with the greatest weight being given to reliabil- 
ity and simplicity of service. After full consideration, it was 
decided to employ the principle of self-refrigeration, where 
carbon dioxide would be discharged in the case of emergency 
or long periods of transit, and the latent heat of vaporization 
required, extracted from the lading through the medium of 
heat exchange. This was achieved by locating two liquid CO» 
containers, approximately 330 Imperial gallons each, inside the 
main tank shell. The use of two containers was decided upon, 
firstly, to supply adequate heat exchange area and, secondly, 
two separate systems to ensure one system always being in 
operation. 

The location of the tanks was arranged so that the liquid 
interface in the CO, tanks and the liquid interface of the lading 
were on similar levels. The internal CO. containers were 
designed to withstand an external pressure in excess of 495 p.s.1., 
the flow rating pressure of the safety valve on the main tank, 
which resulted in considerable over- -design for the internal pres- 
sure. These tanks were hydrostatically tested to 600 p.S.l. prior 
to installation and, of course, are made of the same low temper- 
ature material as the main tank. 


Refrigeration 

The vapor pressure curves of liquid CO, and HCI anhydrous 
bear a close relationship which makes COz an ideal refrigerant 
in this case. The liquid COs» tanks are supplied with pressure 
regulators which are set to discharge at 300 p.s.i. equivalent to 
carbon dioxide at about minus 2°F. By loading the commodity 
at minus 50°F. about 10 days will elapse i in an ambient tempera- 
ture of 90°F. before the discharge temperature of the COsz is 
reached. There is sufficient carbon dioxide available to maintain 
approximately O°F. in this car for an additional 7 days and there 
is adequate time available in the next 20 days before the pressure 
in the tank would approach the safety valve setting during which 
the carbon dioxide containers can be recharged by pumping 
liquid COs from tank truck or tank car at a pressure of about 
275 p.s.i. gauge. 

In the event of leakage of HCI safety valve or angle valves, 
the Type “C” emergency kits av ailable for chlorine tank cars 
can be used without modification to control leakage. (see 
Chlorine Institute Inc., “Types and Location of Chlorine 
Emergency Kits in North America”, First Issue April 1959). 


Construction Details 


Fittings for hydrogen chloride are basically the same as 
fittings used on chlorine cars plus a thermometer well and two 
level telltales. Fi ittings include: 


(1) Crosby safety valve with frangible disc, set at 450 p.s.i. 

(2) Two 1-in. angle valves for liquid control, on the longitudinal 
centre line, with excess flow valves, and 114-in. pipes 
extending to a sump in the bottom of the tank. Pipes are 
screwed to permit maintenance of the excess flow valves, 
as for chlorine cars. 

Two I-in. angle valves for vapor control on the transverse 

centre line. 

Thermometer well to check lading temperature. Closure is 

an approved fitting that permits testing for leaks without 

complete removal. 
(5) Two 3¢-in. angle valves of stainless steel with pipes extend- 
ing into the tank to indicate liquid level at —20°F. and 

— 50°F, 

Fittings for carbon dioxide, pressure regulator and safety 
valve, are the same as those used on carbon dioxide tank cars 
on the theory that the heat flow will be relatively the same as 
a carbon dioxide car and they will be capable of handling the 
same volume. These include the following: 


(1) Safety valve set at 375 p.s.1. 
(2) Pressure regulator set at 300 p.s.i. (equivalent to carbon 
dioxide at —2°F.) 
(3) One I-in. angle valve, same as above, for liquid control, 
with I-in. pipe to bottom of container. 
One I-in. angle valve, same as above, with 1-in. pipe 
extending 134-in. into container to control vapor and also 
act as a fixed length dip tube gauging device for loading of 
carbon dioxide to 95% of volume of container and nozzle. 
Thermometer well to check carbon dioxide temperature. 
‘Tank car shipments were commenced several months ago 
and are now being made on a routine basis in cars owned by 
Canadian General Transit Co. Limited. The temperature well 
and liquid level indicators allow the contents as loaded to be 
checked for level and temperature prior to shipment and in the 
case of incorrect loading, adjusted through the liquid lines. 
Basically simple, the tank car meets the full requirements of 
safety, it is reliable and it provides a container of reasonable 
size and reasonable net weight to ship an economical quantity 
of hydrogen chloride anhy drous by rail. 
While a limited number of cars have been put into service 
in the United States, Dow is the first and so far only Canadian 
shipper of HCI anhydrous by tank car. 
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